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ABSTRACT 


X 


A  new  algorithm  for  error  detection  and  smoothing,  EDS,  has  been  deve¬ 
loped  for  infrared  image  data.  Numerical  validation  indicates  that  EDS  de¬ 
tects  more  than  99%  of  all  significant  errors.  Operating  in  the  time  domain, 
EDS  exploits  smoothness  criteria  calculated  from  laboratory  point  source  meas¬ 
urements,  and  isolates  error  bursts  via  a  region-growing  approach.  EDS  is  ap¬ 
plicable  to  data  from  any  optical  scanner  whose  output  is  sampled  digitally 
multiple  times  per  dwell. 
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1.  INTRODUCTION 


The  Navy’s  Background  Measurements  and  Anal¬ 
ysis  Program  (BMAP)  was  recently  begun  to  obtain 
the  data  needed  to  design  surface-  and  air-based  in¬ 
frared  search  and  track  (IRST)  devices.1,2  The  data 
are  intended  to  be  a  test  set  of  infrared  (IR)  back¬ 
ground  images  for  use  in  off-line  simulation  of  alter¬ 
native  signal  processing  techniques  for  false-alarm 
suppression.  The  IR  backgrounds  data  are  acquired 
with  a  sensor  that  has  high  spatial  resolution,  high  sen¬ 
sitivity,  and  simultaneous  4  to  5  /«m  and  8  to  11  *im 
coverage. 3 

Initial  IR  clutter  measurements  were  obtained  at 
Montauk  Point,  Long  Island,  N.  Y. 1-5  Unfortunate¬ 
ly,  the  initial  measurements  were  contaminated  by  a 
malfunctioning  digital  recorder  with  about  4%  defec¬ 
tive  data  (3%  small-amplitude  errors  and  1%  large- 
amplitude  errors).  Thus,  each  data  frame  of  6400  sam¬ 
ples  contains  (on  the  average,  prior  to  data  reduction) 
about  60  large-amplitude  errors,  i.e.,  errors  of  am¬ 
plitude  10  times  the  sensor  noise  or  greater.  Visual  in¬ 
spection  of  the  data  shows  that  the  errors  tend  to  occur 
in  bursts,  i.e.,  as  multiple  sequential  errors. 

Subsequent  to  the  Montauk  Point  measurements 
there  was  a  hiatus  in  field  activities  during  which  soft¬ 
ware  techniques  were  developed  for  extracting  a  limited 
interim  data  set  for  use  until  better  data  were  ob¬ 
tainable. 

The  need  to  recognize  “wild  points”  or  “outliers” 
in  data  is  a  requirement  of  data  reduction  in  all  ex¬ 
perimental  sciences,  including,  for  example,  chemis¬ 
try,  speech  analysis,  spectroscopy,  and 
geophysics.6'10  The  Montauk  data  analysis  is  unusual 
relative  to  other  outlier  identification  problems  in  that 
the  Montauk  data  contain  a  large  number  of  multi¬ 
ple  sequential  errors.  On  the  average,  each  frame  of 
Montauk  data  contains  one  error  burst  of  10  (or  more) 
sequential  errors,  as  well  as  many  error  artifacts  of 
lesser  length  (Table  I). 


'R.  L.  Lucke,  A.  P.  Schaum,  J.  C.  Kershenstein,  J.  Michalo- 
wicz,  B.  V.  Kessler.  A.  B.  Blumenthal,  I.  Goldstein,  and  A. 
Krutchkoff,  “The  Navy’s  Infrared  Background  Measurement 
and  Analysis  Program,”  Proc.  IRIS  Specially  Group  Meeting 
on  Targets,  Backgrounds,  and  Discrimination,  San  Diego,  11-14 
Feb.  1986. 

2R.  A.  Steinberg,  “Navy  IR  Background  Measurements  and 
Analysis  Program,”  Proc.  Tri-Services  Infrared  Backgrounds 
Symp.,  Oct.  18-20,  1983,  Mitre  Corporation,  pp.  216-223. 

’M.  S.  Longmire,  "A  Final  Technical  Report  on  Calibration  and 
Use  of  Clutter  Data  for  Simulation,”  Western  Kentucky  Uni¬ 
versity,  Bowling  Green,  Ky.,  Contract  No.  N0014-84-C-2034  (un¬ 
dated  draft  report  to  Naval  Research  Laboratory). 


Table  1 

Number  of  extended  artifacts  (i.e.,  multiple  errors)  in 
an  average  frame  of  Montauk  Point  data.  By  “length 
of  artifact”  is  meant  the  number  of  sequential  defec¬ 
tive  samples  in  an  error  burst.  Thus,  an  artifact  of 
length  8  is  an  error  burst  consisting  of  8  sequential 
defective  samples.  Each  data  frame  contains  6400 
samples. 


Length  of  artifact 

(samples),  L  6  7  8-9  s  10 

Average  number  of  1  1  1  1 

artifacts  of  length 
L  (per  frame) 


Identification  of  multiple  sequential  and  nearly  se¬ 
quential  errors  is  more  difficult  than  identification  of 
isolated  single-sample  errors.  Figure  1  illustrates  the 
difficulties  experienced  in  using  a  commonly  available 
outlier  code,  the  IMSL  routine  ICSMOU,  against  a 
pair  of  closely  spaced  single-sample  errors.  In  the  sim¬ 
ple  example  of  Fig.  1 ,  ICSMOU  damages  a  good  data 
point.  Figure  2  illustrates  a  compressed  notation  for 
representing  the  results  of  algorithm  performance  on 
short  sequences  of  synthetically  damaged  data.  The 
single  graph  in  Fig.  2  conveys  the  same  information 
as  the  three  graphs  in  Fig.  1 . 


4A.  Hirschman,  “BMAP  Surface-Based  Background  Measure¬ 
ment  Activities  at  Montauk  Point,  New  York,  August  1983,” 
Naval  Surface  Weapons  Center,  White  Oak,  Silver  Spring,  Md., 
draft  report,  Nov.  1984. 

’D.  R.  Jensen,  “Data  Report  for  IR  BMAP  Experiment,  Aug. 
1983,”  Naval  Ocean  Systems  Center,  San  Diego,  Cal.,  Memoran¬ 
dum  DRJ:plk,  Serial  5325/109,  14  Nov.  1983. 

‘J.  R.  Rice,  Numerical  Methods,  Software,  and  Analysis, 
McGraw-Hill,  New  York,  1983. 

7H.  Ney,  “A  Dynamic  Programming  Technique  for  Nonlinear 
Smoothing”,  Proc.  IEEE  International  Corf.  Acoustics,  Speech 
and  Signal  Processing,  Vol.  I,  Mar.  30-31  and  Apr.  1, 1981,  At¬ 
lanta,  Ga.,  pp.  62-65. 

'H.  Ney,  “A  Dynamic  Programming  Algorithm  for  Nonlinear 
Smoothing,”  Signal  Processing,  Vol.  5,  No.  2,  Mar.  1983,  pp. 
163-173. 

,F.  Pasian  and  A.  Crise,  "Restoration  of  Signals  Degraded  by 
Impulse  Noise  by  Means  of  a  Low-Distortion  Non-Linear  Fil¬ 
ter,"  Signal  Processing,  Vol.  6,  No.  1,  Jan.  1984,  pp.  67-76. 

,0R.  E.  Boucher  and  J.  P.  Noonan,  “Adaptive  Detection  and 
Removal  of  NonGaussian  Spikes  from  Gaussian  Data,”  IEEE 
Trans.  Pattern  Analysis  and  Machine  Intelligence,  Vol.  PAMI-4, 
No.  2,  Mar.  1982,  pp.  132-136. 
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Figure  1 .  Inadequate  treatment  of  multiple  errors  can 
make  the  data  worse,  (a)  Original  data  (constant  val¬ 
ue).  (b)  Damaged  data,  with  outliers  at  samples  4  and  6. 
(c)  “Corrected”  data,  after  application  of  IMSL  routine 
ICSMOU,  designed  for  generic  wild  point  removal  ap¬ 
plications*  ICSMOU  fails  to  identify  data  errors  at 
samples  4  and  6,  and  at  sample  5  damages  what  origi¬ 
nally  had  been  a  good  data  value. 


A  modified  median  filter  approach9  and  nonlinear 
smoothing  by  dynamic  programming7'8  both  appear 
to  hold  some  promise.  However,  both  these  techniques 
require  setting  thresholds  and  selecting  smoothness 
parameters  that  depend  in  a  non-obvious  way  on  the 
character  of  the  data.  The  original  error  detection  tech¬ 
nique  described  in  this  report  is  based  on  exploiting 
smoothness  criteria  that  derive  directly  from  the  known 
point  response  characteristics  of  the  infrared  scanner. 
The  values  ascribed  to  the  smoothness  parameters  are 
simply  calculated  from  measurements  made  against 
a  laboratory  point  source,  and  do  not  depend  on  the 
highly  variable  structure  of  the  1R  clutter  data  viewed 
in  the  field. 

The  computation  is  structured  such  that  the  data 
are  pre-screened  with  a  simple  algorithm  that  allows 
further  attention  to  be  concentrated  only  on  “suspect” 
neighborhoods.  There  are  no  computation-intensive 


Figure  2.  Illustration  of  compressed  notation  for 
representing  the  results  of  algorithm  performance  on 
short  sequences  of  synthetically  damaged  data.  The 
single  graph  in  this  figure  conveys  the  same  informa¬ 
tion  as  the  three  graphs  in  Fig.  1. 0’s  are  original  un¬ 
damaged  data;  vertical  arrows  indicate  data  errors  that 
we  intentionally  implant;  X's  are  values  after  applica¬ 
tion  of  correction  algorithm  (ICSMOU  in  this  instance). 
Where  X  and  0  coincide,  only  O  is  drawn.  Where  X  falls 
at  the  tip  of  an  arrowhead  (e.g.,  samples  4  and  6)  the 
correction  software  has  failed  to  recognize  a  data  er¬ 
ror.  Where  X  appears  at  a  sample  with  no  vertical  ar¬ 
row  (i.e.,  with  no  error),  correction  software  has 
damaged  a  good  data  point  (e.g.,  sample  5).  When  X 
fails  anywhere  other  than  at  the  tail  of  an  arrow,  cor¬ 
rection  software  has  identified  the  data  sample  as  an 
error  but  has  generated  an  imperfect  correction  (not 
illustrated  on  the  figure).  The  dashed  line  is  the  “cur¬ 
rent  waveform,”  i.e.,  the  waveform  as  it  would  appear 
after  the  last  processing  operation  performed  on  the 
data. 


steps,  such  as  FFT’s  or  histograms,  required  for  the 
data  set  as  a  whole.  Execution  time  on  a  VAX  11/780 
was  about  10  times  faster  than  ICSMOU. 

Outlier  identification  and  replacement  is  the  first 
step  in  the  IR  clutter  data  reduction  process.  The  se¬ 
cond  step  in  data  reduction  involves  correction  of  a 
structural  distortion  in  the  longwave  data  caused  by 
resistive  coupling  between  the  longwave  detector  ele¬ 
ments.  "•'2  The  remaining  steps  in  data  reduction  in¬ 
clude  stagger  offset  compensation,  responsivity  nonun¬ 
iformity  compensation,  and  radiometric  calibration.  J 

"R.  A.  Steinberg,  “Infrared  Background  Sensor  Characteriza¬ 
tion  and  Data  Reduction,”  IRIS  Specialty  Group  on  Targets, 
Backgrounds,  and  Discrimination,  Feb.  1-2,  1984,  Camp  Pen¬ 
dleton,  Cal.,  pp.  209-231. 

I2R.  A.  Steinberg,  “Elimination  of  Sensor  Artifacts  from  In¬ 
frared  Data,"  Report  8861 ,  Naval  Research  Laboratory,  Wash¬ 
ington,  D.C.,  11  Dec.  1984. 
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2.  OUTLIERS  APPEAR  IN  THE  DATA 


Prior  to  the  field  tests  at  Montauk  Point,  labora¬ 
tory  bar-target  and  point-target  measurements  were 
obtained  and  analyzed  to  characterize  the  1R  scanner’s 
transient  response  and  radiometric  sensitivity. 11 
Figure  3  shows  a  sample  of  this  data,  obtained  against 
a  laboratory  target  consisting  of  six  vertical  bars  of 
progressively  increasing  width. 

The  scan  line  shown  in  Fig.  3  contains  three  single¬ 
sample  errors,  a  large  positive-going  error  at  sample 
315,  a  very  large  positive-going  error  at  sample  422, 
and  a  small  negative-going  error  at  sample  433.  Sen¬ 
sor  closed-cover  noise  is  about  one  count,  RMS;  major 
y-axis  divisions  are  500  counts,  as  shown  on  the  figure. 
Thus,  none  of  the  waveform  structure  visible  in  Fig. 
3  is  attributable  to  sensor  closed-cover  noise. 

The  narrowest  bar,  left-most  in  Fig.  3,  is  spatially 
unresolved.  The  physically  smallest  possible  optical 
source,  viz.,  an  ideal  point  source,  induces  a  pulsed- 
shaped  waveform  (“the  point  response”)  nearly  iden¬ 
tical  to  this  narrowest  bar.  Any  pulsed-shaped  wave¬ 
form  narrower  than  the  point  response  must  be  an 
artifact  of  electrical,  rather  than  optical,  origin.  Thus, 
applying  a  simple  pulse-width  criterion,  single-sample 
errors  are  detected  easily  both  by  eye  and  by  algorithm. 
Intuitively,  the  key  to  simple  identification  of  single¬ 
sample  errors  lies  in  prior  knowledge  of  the  system 
point  response. 


Figure  4  depicts  a  scan  line  of  data  from  Montauk 
Point,  obtained  looking  up  from  the  surface  against 
a  cloud/sky  background.  The  data  shown  in  Fig.  4 
contain  about  three  times  as  many  multiple  errors  as 
an  average  scan  line:  of  the  400  data  samples  that  com¬ 
prise  the  waveform,  24  are  multiple  errors.  Three  badly 
damaged  data  strings,  each  coincidentally  of  1 1  sam¬ 
ples  duration,  are  found  to  contain  26  of  the  30  total 
errors.  The  first  badly  damaged  11-sample  interval, 
centered  on  sample  171,  contains  7  errors  and  4  good 
data  points  (Fig.  5a).  The  second  error  burst,  centered 
on  sample  295,  contains  8  errors  and  3  good  data 
points  (Fig.  5b).  The  third  and  final  error  burst,  cen¬ 
tered  on  sample  422,  contains  1 1  errors  and  no  good 
data  points  (Fig.  5c). 

The  remaining  4  errors  in  the  scan  line  are  each  iso¬ 
lated  single-sample  errors.  These  isolated  outliers,  once 
identified  (detected),  can  be  replaced  by  interpolation 
without  significantly  degrading  data  quality,  due  to 
the  high  data  sampling  rate  (3.47  samples  per  dwell). 

Although  the  IR  clutter  data  within  the  3  multi-error 
regions  cannot  be  reconstructed,  this  only  marginally 
affects  the  utility  of  the  data.  Once  the  3  error  bursts 
are  identified  and  excised,  the  remaining  data  may  then 
be  used  for  a  variety  of  IRST  processor  simulation 
studies. 

We  note  that  our  notation  on  “zoom”  plots  of  er¬ 
rors  found  in  field  test  data  (e.g..  Fig.  5)  is  different 
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Figure  3.  Single  scan  line  (midwave  channel  15)  ob-  Figure  4.  Scan  line  of  badly  damaged  data.  The  data 

tained  against  a  laboratory  target  consisting  of  six  ver-  shown  in  this  figure  contain  about  four  times  as  many 

tical  bars  of  progressively  increasing  width.  Single-  errors  as  an  average  scan  line:  of  the  400  data  sam- 

sample  errors  occur  at  sample  numbers  315,  422,  and  pies  that  comprise  the  waveform,  24  are  multiple  er- 

433.  A  vertical  tick  mark  above  the  abscissa  indicates  rors.  Vertical  lines  connect  data  errors  to  interpolated 

the  location  of  each  data  error.  values  on  the  waveform. 
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Figure  5.  Expanded  view  of  the  three  error  burst  intervals  in  Fig.  4.  Circles  not  connect¬ 
ed  to  vertical  lines  indicate  valid/undamaged  data  points.  Vertical  lines  connect  interpo¬ 
lated  values  (indicated  by  circles)  to  damaged  data  points.  Arrowhead  indicates  that  error 
value  is  off-scale,  (a)  11 -sample  interval  centered  on  sample  171  contains  7  errors,  (b) 
11 -sample  interval  centered  on  sample  295  contains  8  errors,  (c)  11-sample  interval  cen¬ 
tered  on  sample  422  contains  11  errors. 
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than  for  data  into  which  we  ourselves  have  implant¬ 
ed  errors  for  algorithm  evaluation  purposes  (“synthetic 
errors”,  e.g.,  Fig.  2).  Our  notation  for  error  detec¬ 
tion  in  field  test  data,  illustrated  in  Fig.  5,  is  as  follows: 

•  Vertical  lines  indicate  samples  identified  as  errors 
by  the  EDS  software. 

•  Open  circles  attached  to  vertical  lines  are  inter¬ 


polated  values,  i.e.,  replacement  values  for  iden¬ 
tified  errors. 

•  Unattached  ends  of  vertical  lines  are  original, 
presumably  defective,  values. 

•  Arrow  heads  indicate  off-scale  error  values. 

•  Open  circles  not  attached  to  vertical  lines  indi¬ 
cate  “good”  measured  data. 
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3.  HOW  RELIABLY  ARE  THE  ERRORS  DETECTED? 


The  example  of  Figs.  4  and  5  illustrates  the  capa¬ 
bility  of  our  new  error  detection  algorithm  for  isolat¬ 
ing  long  sequences  of  damaged  data.  In  this  section 
we  present  the  results  of  a  validation  analysis  intended 
to  address  the  following  statistical  aspects  of  algorithm 
performance: 

•  What  fraction  of  data  errors  elude  detection  by 
EDS?  How  many  of  the  undetected  errors  are  sig¬ 
nificant,  in  the  sense  of  having  large  amplitude? 

•  What  fraction  of  the  good  data  is  damaged  by 
EDS? 

•  For  those  data  errors  that  are  correctly  identified 
as  errors  by  EDS,  how  accurate  are  the  applied 
corrections? 

Our  approach  to  statistical  validation  was  first 
described  in  broad  outline  in  Ref.  2. 

The  validation  process  starts  with  the  selection  of 
a  data  set  A  that  is  error-free  by  the  standard  of  the 
software  error  detection  and  smoothing  (EDS)  oper¬ 
ator,  0,  represented  as 

A  =  0(A),  (1) 

i.e.,  operating  on  data  set  A  with  the  EDS  software 
leaves  the  data  set  unchanged. 

The  data  set  A  used  for  analysis  consisted  of  800 
scan  lines,  such  as  that  plotted  in  Fig.  6. 

We  numerically  generate  and  add  to  A  a  set,  S,  of 
error  values,  in  which  every  nonzero  component 


Sample  number 

Figure  8.  An  error-free  scan  line.  800  such  error-free 
lines  of  data  (collectively  referred  to  as  data  set  A)  are 
used  in  connection  with  the  statistical  validation  anal¬ 
ysis  of  the  error  detection  algorithm.  The  significance 
of  regions  labeled  A  to  E  is  discussed  in  the  text. 


represents  a  synthetic  data  error  of  known  location 
and  amplitude  (Fig.  7).  The  result  is  a  synthetically 
corrupted  data  set  B  (Fig.  8), 

B  =  A  +  S  .  (2) 
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Figure  7.  Sample  scan  line  of  synthetic  data  errors, 
containing  two  five-fold  errors  (A  and  C),  two  double 
errors  (B  and  D),  and  an  isolated  error  (E).  The  error 
generation  algorithm  was  designed  to  provide  error 
statistics  specified  with  respect  to  an  amplitude  pro¬ 
file  (Fig.  9)  and  the  number  of  error  bursts  of  a  given 
duration  (Fig.  10).  The  ensemble  of  800  such  error  lines 
is  referred  to  in  the  text  as  data  set  S. 


Sample  number 

Figure  8.  A  synthetically  corrupted  scan  line  obtained 
by  adding  synthetic  data  errors  (Fig.  7)  to  an  error-free 
scan  line  (Fig.  6).  Error  bursts,  i.e.,  multiple  sequen¬ 
tial  errors,  are  injected  at  positions  denoted  as  A,  B, 
C,  and  D,  on  the  waveform;  a  single  error  was  injected 
at  position  E.  The  multiple  error  at  D  occurs  on  a  peak, 
as  seen  by  comparison  with  the  corresponding  region 
D  in  Fig.  6. 
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The  error  set,  S,  is  described  by  two  histograms 
(Figs.  9  and  10)  that  statistically  characterize  the  syn¬ 
thetic  errors  with  respect  to  amplitude  and  spatial  ex¬ 
tent,  and  that,  taken  together,  constitute  a  model  for 
the  error  generation  process. 

Each  data  frame  used  in  the  statistical  validation 
had  6400  samples: 

16  scan  lines  per  frame  x  400  samples  per  line 
=  6400  samples  per  frame. 

Fifty  statistically  independent  error  frames  were  nu¬ 
merically  generated,  containing  a  total  of  1 1 ,470  data 
errors.  It  follows  that  the  overall  error  rate  is: 

(11,470  defective  samples) 
error  rstc  =  — . . .  —  -  — . 

6400  samples  per  frame  x  50  frames 

=  0.0358  =  3.58%. 

Figure  9  shows  how  the  1 1 ,470  synthetic  errors  were 
distributed  in  amplitude.  For  example,  8458  errors 
(about  three-fourths  of  the  total)  had  amplitude  values 
between  1  and  9  counts.  The  number  of  errors  of  am¬ 
plitude  >  10,  as  a  fraction  of  the  total  data  set,  is  thus 
given  by: 

((11,470  -  8,458)/320,000]  =  0.94%. 
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Figure  9.  Amplitude  histogram  of  1 1 ,470  synthetical  ly 
generated  errors,  dispersed  into  50  frames  of  data 
(overall  error  rate  =  3.58%).  About  three-fourths  of  the 
errors  have  values  from  1  to  9  counts  in  amplitude.  The 
data  error  sign  was  modeled  as  statistically  Indepen¬ 
dent  of  the  error  amplitude,  and  was  assigned  positive 
and  negative  values  with  equal  likelihood.  The  synthet¬ 
ic  errors  are  generated  by  a  program  that  closely 
matches  their  amplitude  histogram,  shown  here,  to  the 
amplitude  histogram  of  errors  found  in  Montauk  Point 
data. 
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Figure  10.  Artifact  length  histogram  of  synthetically 
generated  errors.  A  total  of  7166  error  artifacts  (i.e.,  er¬ 
ror  bursts)  were  dispersed  into  50  data  frames.  As  in¬ 
dicated  by  the  left-most  bar  in  the  figure,  about  82% 
of  the  artifacts  were  single-sample  errors.  The  synthet¬ 
ic  errors  are  generated  by  a  program  that  closely 
matches  their  artifact  length  histogram,  shown  here, 
to  the  artifact  length  histogram  of  errors  found  in  Mon¬ 
tauk  Point  data.  Figures  9  and  10  together  provide  a 
characterization  of  the  Montauk  data  as  well  as  a 
model  for  the  synthetic  error  generation  process. 


We  consequently  speak  of  the  data  as  having 
suffered  about  1%  significant  damage. 

In  this  report  a  group  of  contiguous  errors  bound¬ 
ed  on  each  side  by  a  good  data  value  is  referred  to 
either  as  an  error  burst  or  as  an  error  “artifact,”  with 
the  artifact  length  being  defined  as  the  number  of  se¬ 
quential  errors.  The  1 1,470  errors  in  error  set  S  were 
grouped  into  a  total  of  7166  artifacts,  with  the  distri¬ 
bution  of  lengths  as  shown  in  Fig.  10.  For  example, 
Fig.  10  indicates  that  5864  artifacts  (82%  of  the  total) 
were  of  unit  length,  i.e.,  were  isolated  single-sample 
errors. 

When  our  EDS  algorithm  is  applied  to  the  synthetic¬ 
ally  damaged  waveform  shown  in  Fig.  8,  the  smoothed 
waveform  shown  in  Fig.  1 1  is  obtained.  Comparing 
visually  the  original  and  smoothed  waveforms  (Figs. 
6  and  1 1 ,  respectively)  we  note  that  all  regions  of  data 
damage  have  been  identified  and,  with  varying  degrees 
of  success,  corrected.  The  multiple  errors  at  B  and  C 
and  the  single  error  at  E  have  been  corrected  with  high 
accuracy.  The  multiple  error  at  A  has  been  imperfectly 
corrected,  with  the  result  that  the  trough  present  at 
A  in  the  original  data  (Fig.  6)  has  been  sharpened 
somewhat  in  the  corrected  data  (Fig.  II).  The  largest 
residual  effect  of  the  damage/smoothing  process  ap¬ 
pears  at  location  D.  The  peak  originally  present  at  D 
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Figure  1 1 .  Smoothed  scan  line  obtained  by  applying 
Algorithm  EDS  to  the  synthetically  damaged  scan  line 
depicted  in  Fig.  8.  Neighborhoods  labeled  A  to  E  are 
compared  in  the  text  of  the  report  with  corresponding 
neighborhoods  in  Fig.  6. 


in  Fig.  6,  badly  damaged  by  a  multiple  error  (Fig.  8), 
has  been  sliced  away  by  EDS  (Fig.  1 1).  The  data  er¬ 
rors  at  D  were  successfully  detected  but  badly  cor¬ 
rected. 

Since  restoration  of  the  true  data  structure  is  gener¬ 
ally  not  possible  within  regions  where  multiple  errors 
have  occurred,  the  EDS  algorithm  is  better  described 
as  providing  error  detection  and  censoring,  rather  than 
error  correction,  within  the  “blind  spots”  created  by 
multiple  errors. 

Considering  the  example  of  Figs.  6,  8,  and  1 1,  the 
net  effect  of  data  damage/smoothing  has  been  that 
one  structural  feature  in  the  original  data  has  been  lost 
(D),  another  local  feature  has  been  somewhat  modi¬ 
fied  in  structure  (A),  and  all  other  features  have  been 
either  repaired  or  left  unchanged  in  their  original  un¬ 
damaged  condition.  Whether  the  smoothed  data,  Fig. 
1 1 ,  provides  a  sufficiently  accurate  representation  of 
the  original  data.  Fig.  6,  is  apparently  application- 
dependent. 

We  next  consider  error  correction  statistics  estab¬ 
lished  over  800  scan  lines  (i.e.,  50  frames)  of  synthet¬ 
ically  damaged  data.  In  the  terminology  of  Eqs.  1  and 
2,  we  operate  on  the  synthetically  corrupted  data  set, 
B,  with  algorithm  0, 

A  =  e(B) ,  (3) 

to  obtain  an  estimate,  A,  for  the  original  error-free 
data  set,  A.  In  our  numerical  study,  an  ensemble  of 
800 corrected  waveforms  (e.g..  Fig.  1 1)  constitutes  A. 

Our  numerical  study  uses  data  sets  A,  S,  and  A  that 


each  consist  of  50  frames,  i.e.,  320,000  values.  Equiva¬ 
lently,  we  may  consider  corresponding  elements  in 
these  data  sets  as  forming  320,000  3-tupies,  ( a.s.a ), 
where  a  denotes  a  single  element  of  A  (and  similarly 
for  s  and  a).  A  test  is  performed  on  each  3-tuple,  with 
four  possible  outcomes,  as  indicated  in  Table  2. 

Considering  an  Outcome  2  result,  for  example,  we 
see  from  Table  2  that  this  corresponds  to  the  follow¬ 
ing  conditions  on  the  3-tuple  ( a.s.a ): 


s  =  0 

(4a) 

a  *  a 

(4b) 

Equation  4a  indicates  that  the  data  sample  in  ques¬ 
tion  has  been  modified  with  a  zero-amplitude  error, 
i.e.,  the  data  sample  is  a  valid  data  point.  Equation 
4b  indicates  that  the  processed  data  value,  3,  is  not 
equal  to  the  original  data  value,  a.  The  interpretation 
of  the  Outcome  2  result  is,  thus,  straightforward:  a 
good  data  value  has  been  incorrectly  identified  as  an 
error  by  the  EDS  algorithm  and  modified  in  ampli¬ 
tude,  i.e.,  “damaged,”  by  the  amount  3  -  a. 

Of  the  320,000  trials,  an  Outcome  2  result  is  ob¬ 
tained  618  times,  i.e.,  618  valid  data  points  are 
damaged.  The  618  corresponding  values  of  13  -  a  I 
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Figure  12.  Amplitude  histogram  of  damaged  points. 
Of  the  320,000  sample  values  (i.e.,  50  frames)  in  the  vali¬ 
dation  data  set  A,  618  valid  data  points  are  incorrectly 
classified  as  errors  and  are  subsequently  modified  in 
amplitude  ("damaged”)  by  the  EDS  Algorithm.  The  left¬ 
most  bar  of  the  histogram  indicates  that  480  of  the 
damaged  samples  (78%)  were  damaged  by  4  or  less 
digital  counts.  Good  data  values  damaged  by  more 
than  19  counts  were  located  at  the  edge  of  an  error 
burst. 
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are  depicted  as  a  histogram  in  Fig.  12.  If,  for  a  given 
application,  data  damage  of  less  than  10  counts  is 
deemed  insignificant,  it  follows  from  Fig.  12  that  27 
data  points  have  been  “significantly  damaged”  by  the 
EDS  algorithm. 


From  Eqs.  5a  and  6  it  follows  that  Outcome  1  cor¬ 
responds  to 

JKO  (7a) 

a  =  b  .  (7b) 


Table  2 

A  test  is  performed  on  each  3-tuple,  ( a,s,i )  with  four 
possible  outcomes.  Since  each  data  set,  (A,  S,  and 
A)contains  320,000  values  (corresponding  to  50 
frames),  this  test  was  performed  320,000  times.  The 
values  of  Is  I  corresponding  to  outcome  1  (missed  er¬ 
rors)  were  used  to  assemble  the  histogram  in  Fig.  13. 
Similarly,  values  of  Id  -  a  I  corresponding  to  out¬ 
comes  2  and  3  were  used  to  assemble  Figs.  12  and  14, 
respectively. 


Number  of 

Outcome  s  a  -  a  Interpretation  Occurrences 

1  * 0  5  Bad  data  has  2,839 

not  been 
detected 

2  0*0  Good  data  has  618 

been  damaged 

3  *0  *s  Data  error  has  8,631 

been  detected 

4  0  0  Good  data  has  307,912 

been  left 
unchanged 


The  interpretation  of  Table  2  for  Outcome  1  is  only 
slightly  less  obvious  than  Outcome  2.  According  to 
Table  2,  Outcome  1  corresponds  to 

5*0  (5a) 

&  -  a  =  s  (5b) 


From  Eqs.  5b  and  2, 


<?  =  a  +  s  =  b  ,  (6) 

where  b  is  a  sample  from  the  damaged  data  set,  B. 


Equation  7a  indicates  that  the  data  sample  in  ques¬ 
tion  has  been  modified  by  a  nonzero  value,  i.e.,  the 
data  sample  is  an  error.  Equation  7b  indicates  that 
the  processed  data  value,  a,  is  exactly  equal  to  the 
defective  data  value  b,  i.e.,  no  correction  has  been  per¬ 
formed.  Thus,  our  interpretation  of  an  Outcome  1 
result  is  that  a  data  error  of  amplitude  s  has  propagated 
unchanged  through  EDS,  i.e.,  the  error  has  been 
“missed.” 

Of  the  320,000  trials,  an  Outcome  1  result  is  ob¬ 
tained  2839  times,  i.e.,  2839  errors  are  missed.  The 
2839  corresponding  values  of  Is  I  are  depicted  as  a 
histogram  in  Fig.  13.  If,  for  a  particular  application, 
data  errors  of  less  than  10  counts  are  deemed  insig¬ 
nificant,  it  follows  from  Fig.  13  that  24  “significant¬ 
ly  large”  errors  have  not  been  detected  by  EDS. 
However,  we  see  from  Fig.  9  that  3012  errors  were 
implanted  with  amplitudes  a  10.  Thus,  the  fraction 
of  significantly  large  errors  not  detected  by  EDS  is 
(24/3012)  =  0.8%;  slightly  better  than  99  out  of  each 
100  large  amplitude  errors  are  detected. 
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Figure  1 3.  Amplitude  histogram  of  undetected  errors. 
Of  the  1 1,470  errors  Implanted  into  50  frames  of  data, 
2839  errors  are  undetected  by  EDS.  However,  only  24 
of  the  missed  errors  were  of  significant  amplitude,  i.e, 
amplitude  10  or  larger.  Thus,  after  correction  there 
were,  on  average,  0.5  significant  undetected  errors  per 
frame.  The  largest-amplitude  undetected  error  (ampli¬ 
tude  39)  occurred  at  the  edge  of  an  error  burst. 
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Similarly,  it  may  be  shown  that  a  Table  2/Outcome 
3  result  corresponds  to  a  trial  in  which  a  data  error 
has  been  correctly  identified  as  an  error  by  EDS,  and 
a  correction  has  been  applied  to  the  data.  Of  the 
320,000  trials.  Outcome  3  occurs  8631  times.  The  8631 
corresponding  values  of  Id  -  a  I  are  represented  as 
a  histogram  in  Fig.  14,  indicating  the  accuracy  of  data 
correction  for  those  data  errors  correctly  identified  as 
errors  by  EDS.  For  example,  if  correction  errors  of 
less  than  20  counts  are  deemed  insignificant,  it  fol¬ 
lows  from  Fig.  14  that  only  21  error  corrections  by 
EDS  were  “significantly  inaccurate.”  The  poorly  cor¬ 
rected  samples  —  the  tail  of  the  Fig.  14  histogram  — 
correspond  to  “sliced  peaks,”  such  as  occurred  at  lo¬ 
cation  D  in  Figs.  6  and  11. 

The  three  histograms  presented  here  as  Figs.  12 
through  14  together  constitute  the  results  of  the  statisti¬ 
cal  assessment/validation  of  our  new  error  detection 
and  smoothing  algorithm,  EDS. 

Additional  comments  regarding  final  assessment  of 
EDS  are  provided  in  the  Conclusions  section. 


Figure  14.  Amplitude  histogram  of  residual  errors  af¬ 
ter  data  smoothing.  Of  the  8631  errors  that  were  cor¬ 
rectly  identified  as  errors  by  EDS,  8078  corrections 
(93.6%)  were  accurate  to  within  4  digital  counts.  Sam¬ 
ples  accumulated  into  the  tail  of  the  histogram  are  at¬ 
tributable  to  “sliced  peaks,"  such  as  occurred  at 
location  D  in  Figs.  6  and  11,  and  multiple  errors  badly 
corrected  at  the  ends  of  a  scan  line. 


4.  SYSTEM  POINT  RESPONSE 


As  discussed  above,  inspection  of  Fig.  3  supports 
the  intuitive  assertion  that  prior  knowledge  of  system 
point  response  allows  us  to  reliably  visually  detect 
single-sample  errors:  since  single-sample  errors  give 
rise  to  spikes  of  width  narrower  than  the  point 
response,  a  pulse-width  criterion  serves  to  detect  single¬ 
sample  errors.  This  intuitive  observation  subsequent¬ 
ly  serves  as  the  basis  for  developing  an  error-detection 
system  for  multiple  errors  as  well  as  single  errors.  Since 
the  key  to  our  approach  lies  in  exploiting  prior  know¬ 
ledge  of  system  point  response  (SPR)  we  will  first  dis¬ 
cuss  the  SPR. 

Each  of  the  IR  scanner’s  32  detectors  has  a  differ¬ 
ent  point  response.  When  we  speak  of  system  point 
response  we  are  referring  to  the  ensemble  of  measured 
point  responses  in  which  each  detector  is  character¬ 
ized  individually. 

The  point  response  measurements  presented  here 
were  obtained  by  horizontally  scanning  the  sensor’s 
field  of  view  across  a  laboratory  point  source.  (The 
laboratory  set-up  is  described  in  Ref.  3.)  The  data  are 


recorded  digitally  with  a  sample  depth  of  12  bits.  Ex¬ 
amples  of  measured  point  responses  are  given  in  Fig. 
IS  for  longwave  channels  S  and  2. 

Inspection  of  Fig.  IS  shows  that  the  point  response 
full-width  at  10%  peak  amplitude  is  about  6.4  sam¬ 
pling  intervals.  The  zero-amplitude  baseline  is  relia¬ 
bly  captured  in  about  12  sampling  intervals,  i.e.,  by 
representing  the  point  response  as  a  13-sample  wave¬ 
form. 

The  sampling  interval  is  96  ^rad  (about  0.55%  of 
1  °)  in  azimuth,  established  by  slaving  the  A/D  con¬ 
verter  to  a  16-bit  shaft  angle  encoder  on  the  sensor’s 
scanning  mirror: 

(2t  rad/rev)/(216  samples/rev)  =  %  prad/sample 

(8) 

The  instantaneous  field  of  view  of  each  detector  is  'A 
mrad  square,  as  determined  by  the  cell  size  (2  mils) 
and  optics  focal  length  (6  inches).  Thus,  the  analog 
detector/postamplifier  output  voltages  are  sampled 
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electronically  at  a  rate  of 

(333  #trad/dwell)/(96  ^rad/sample) 

=  3.47  sampies/dwell.  (9) 

The  shape  of  the  point  response  is  determined  by 
a  number  of  factors  in  addition  to  sample  rate,  in¬ 
cluding 

•  Optics  blur, 

•  Detector  size,  shape,  and  spatial  responsivity 
characteristics, 

•  Postamplifier  frequency  response, 

•  Point  source  position  cross-scan,  relative  to  the 
detector’s  geometric  center,  and 

•  Point  source  position  in-scan  relative  to  the  cell 
geometric  center  and  time-phasing  of  the  A/D 
sampling  clock. 


Figure  IS.  Measured  point  responses  of  longwave 
channel  5  (a)  and  longwave  channel  2  (b).  The  wave¬ 
form  peak  in  (a)  is  very  nearly  an  isosceles  triangle 
while  that  in  (b)  is  very  nearly  flat.  The  difference  in 
structure  between  the  peaks  of  these  two  waveforms 
is  attributable  to  a  small  displacement  of  the  optical 
point  source  in  the  direction  of  scan,  as  discussed  in 
connection  with  Figs.  16  and  17. 


We  note  that  each  detector  is  followed  by  a  wide¬ 
band  preamplifier,  followed  in  turn  by  an  analog 
postamplifier.  The  postampli tiers  each  have  a  1-pole 
iowpass  characteristic  with  a  3  dB  cut-off  at  1000  Hz. 
If  the  postamp  response  were  the  only  pulse  broaden¬ 
ing  mechanism,  the  point  response  in-scan,  p(0),  would 
thus  be  a  simple  exponential, 

fo  .  e < o 

P(0)  =  j 

(^exp (-*„«),  6a0,(10) 

where, 

0  =  azimuth  angle,  samples 
kco  =  cut-off  spatial  frequency,  radians/sample, 

and  the  optical  point  source  is  presumably  located  at 
0  =  0.  The  cutoff  spatial  frequency  is  given  by 

^  -  Itfco/Q  (11) 

where 

fco  -  postamp  3  dB  cutoff  =  1000  sec  ' 1 

0  =  angular  rate  of  scan,  samples/sec. 

The  angular  rate  of  scan  is  calculated  in  appropriate 
measure  (samples/sec)  as: 

6  =  (30.7  degrees/secK17.45  mrad/degree) 
x  (1  sampIe/0.096  mrad)  =  5580  samples/sec. 

(12) 

From  Eqs.  II  and  12, 

kco  =  1.126  radians/sample.  (13) 

If  the  postamp  alone  were  responsible  for  pulse 
broadening,  the  point  response  full-width  at  10%  peak 
amplitude  (denoted  as  0^)  would  be  found  from  Eq. 
10, 

PiO'itc)  =  exp(  )  =  0.1.  (14) 

From  Eqs.  13  and  14, 

0elec  =  2.1 /k co  *  2  samples  (15) 

We  see  from  Fig.  15  that  the  actual  full-width  at  10% 
peak  amplitude,  0fw ,  is  about  6.4  samples.  Assuming 
that  the  separate  contributors  to  the  point  response 
each  contribute  approximately  as  the  square  root  of 
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the  sum  of  squares,  we  find  that  the  postamp  response 
contributes  only  marginally  (about  0.3  samples)  to  the 
overall  width  of  the  point  response: 

0/„.  -  (0/„  -  ohec ) ,  ;  *  0.3  samples. 

We  conclude  from  the  analysis  of  Eqs.  10  to  1 5  that, 
as  expected,  the  measurement  system’s  electrical  fre¬ 
quency  response  is  the  smallest  of  the  three  major  ef¬ 
fects  that,  together,  determine  the  width  of  the  point 
response  (cf.  Table  3). 


Table  3 

Point  response  width  is  determined  by  the  optics  blur, 
the  detector  size,  and  the  postamplifier  frequency 
response.  Tabulated  here  are  the  full-widths  at  10% 
maximum  amplitude  of  the  partial  responses  due  to 
each  of  the  three  main  pulse-broadening  mechanisms. 
Entries  2  and  3  are  from  Eqs.  9  and  15,  respectively. 
Entry  4  is  measured  from  Fig.  15.  Entry  1  is  inferred 
by  assuming  that  the  partial  contributions  add  as  the 
root  of  the  sum  of  squares. 


Full  width  at 
10%  maximum 

_ amplitude 

1.  Optics  blur  (longwave)  5.0  samples 

2.  Detector  size  3.5  samples 

3.  Electrical  frequency  response  2.0  samples 

4.  Total  (longwave  response)  6.4  samples 


Visual  inspection  of  the  point  response  waveforms 
shows  that,  while  their  widths  at  10%  maximum  am¬ 
plitude  are  nearly  the  same,  their  detailed  structure 
near  the  waveform  peak  exhibits  considerable  varia¬ 
tion.  For  example,  we  see  that  the  waveform  peak  in 
Fig.  15a  is  very  nearly  an  isosceles  triangle,  while  that 
in  Fig.  15b  is  very  nearly  fiat.  We  next  show  that  var¬ 
iation  in  the  detailed  structure  near  the  waveform  peak 
is  attributable  to  variation  in  the  in-scan  position  of 
the  optical  point  source. 

An  isosceles  peak  (e.g.,  Fig.  15a)  is  produced  by  the 
sampling  circumstance  depicted  in  Fig.  16.  At  the  bot¬ 
tom  of  Fig.  16  we  schematize  a  circular  optics  blur 
moving  at  a  uniform  speed  across  a  stationary  square 
detector.  The  origin  of  the  time  axis  in  Fig.  16,  t  = 
0,  is  defined  as  the  instant  when  the  optics  blur  is  per¬ 
fectly  centered  on  the  detector,  i.e.,  the  time  when  the 
detected  signal  is  maximum.  Integral  values  on  the 
"sample”  axis  in  Fig.  16  correspond  to  times  when 
the  analog-to-digital  (A/D)  converter  samples  the 


postamp  output.  Figures  1 5a  and  16  correspond  to  the 
unusual  circumstance  where  one  of  the  sampling  times 
occurs  by  coincidence  simultaneous  with  the  reference 
time  /  =  0.  Thus,  a  sample  is  obtained  (sample  6  in 
this  case)  just  as  the  postamp  output  reaches  its  max¬ 
imum  value. 

A  fiat-topped  peak  (e.g..  Fig.  15b)  is  produced  by 
the  sampling  circumstance  depicted  in  Fig.  17.  In  this 
case  the  optics  blur  is  by  coincidence  equally  off-center 
(in  opposite  directions)  at  samples  6  and  7. 


4  5  6  7  8  Sample  number 


Figure  16.  Sampling  conditions  that  result  in  a  sym¬ 
metrically  triangular  (isosceles)  point  response  peak, 
such  as  Fig.  15a.  Quantity  r,  is  the  sampling  interval 
(0.179  msec);  quantity  t„  is  the  dwell  time  (0.622 
msec). 


Most  generally,  the  peak  of  the  point  response  has 
an  asymmetric  triangular  shape,  having  neither  the 
isosceles  triangular  shape  of  Fig.  1 5a  nor  the  fiat  top 
of  Fig.  15b. 

In  performing  the  point  response  measurements 
presented  in  this  section  cross-scan  position  was  ad¬ 
justed  simply  to  maximize  the  point  response  peak  am¬ 
plitude.  No  measurements  were  performed  suitable  for 
assessing  structural  variations  of  the  point  response 
with  sub-pixel  deviations  in  point  source  cross-scan 
position. 
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Figure  17.  Sampling  conditions  that  result  in  a  flat- 
topped  point  response,  such  as  Fig.  15b. 


Figure  19.  Measured  point  responses  of  (a)  midwave 
channel  5  and  (b)  midwave  channel  2. 


Figure  IB.  Average  of  measured  point  responses  for 
(a)  the  longwave  array  and  (b)  the  midwave  array.  Of 
the  16  longwave  detectors,  data  from  only  13  were 
available  for  forming  the  longwave  average  response. 
Data  were  available  from  only  12  of  the  16  midwave 
cells. 


Figure  18  depicts  two  sample  midwave  point 
responses.  The  midwave  point  responses  are  general¬ 
ly  very  similar  to  one  another  and  to  the  longwave 
point  responses. 

Table  4  includes  the  full  complement  of  point 
response  data  available  at  the  time  the  error  detection 
research  was  performed.  The  similarity  in  structure 
of  the  measured  point  responses  in  Table  4  suggests 
that  for  some  applications  it  may  suffice  to  use  an 
array-average  point  response  (Fig.  19).  However,  the 
array  average  responses  are  not  used  or  referred  to 
at  any  later  point  in  this  report. 

The  peak  signal-to-noise  ratios  (S/N)  of  the  point 
response  measurements  provided  in  Table  4  are: 


S/N  * 


f  750,  midwave 
250,  longwave. 


(16) 


f; 
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Table  4 

Normalized  point  responses  for  the  midwave  and  longwave  arrays  used  at  Montauk  Point.  No  point  responses 
and  no  clutter  data  are  available  for  optically  inactive  elements:  midwave  channels  14  and  16,  and  longwave 
channel  1.  Point  response  data  were  unavailable  for  longwave  channel  16.  Point  response  measurements 
were  damaged  by  data  recorder  errors  and  hence  are  not  tabulated  for  midwave  channels  7  and  9,  and  long¬ 
wave  channel  4. 

a.  Normalized  MWIR  Point  Responses 

Channel  1  2  3  4  5  6  7  8  9  10  II  12  13 

1  0.000  0.001  0.073  0.449  0.815  1.000  0.966  0.577  0.244  0.080  0.026  0.008  0.003 

2  0.000  0.005  0.075  0.412  0.786  1.000  0.976  0.610  0.262  0.088  0.030  0.013  0.006 

3  0.000  0.007  0.102  0.467  0.814  1.000  0.941  0.567  0.245  0.085  0.028  0.010  0.004 

4  0.000  0.012  0.171  0.553  0.856  1.000  0.925  0.551  0.233  0.079  0.026  0.011  0.003 

5  0.000  0.012  0.169  0.569  0.873  1.000  0.856  0.458  0.178  0.061  0.021  0.009  0.005 

6  0.000  0.029  0.258  0.641  0.910  1.000  0.824  0.435  0.166  0.054  0.017  0.005  0.003 

8  0.000  0.028  0.244  0.620  0.888  1.000  0.852  0.474  0.193  0.068  0.023  0.008  0.003 

10  0.000  0.027  0.269  0.652  0.913  1.000  0.785  0.400  0.153  0.049  0.016  0.007  0.005 

11  0.000  0.017  0.219  0.603  0.892  1.000  0.796  0.422  0.176  0.057  0.021  0.008  0.003 

12  0.000  0.035  0.307  0.694  0.938  1.000  0.741  0.359  0.131  0.042  0.012  0.005  0.001 

13  0.000  0.030  0.284  0.679  0.937  1.000  0.724  0.355  0.129  0.045  0.015  0.007  0.003 

15  0.000  0.054  0.376  0.762  0.978  1.000  0.669  0.302  0.107  0.033  0.011  0.005  0.003 

Avg  0.000  0.021  0.212  0.592  0.883  1.000  0.838  0.459  0.185  0.062  0.020  0.008  0.003 
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b.  Normalized  IWIR  Point  Responses 

Channel  1  2  3  4  5  6  7  8  9  10  11  12  13 

2  0.000  0.004  0.008  0.051  0.262  0.688  1.000  0.992  0.489  0.165  0.068  0.038  0.013 


3  0.004  0.016  0.033  0.082  0.392  0.890 

5  0.004  0.020  0.035  0.074  0.336  0.836 

6  0.000  0.008  0.051  0.210  0.611  0.946 

7  0.000  0.004  0.019  0.038  0.134  0.622 

8  0.000  0.012  0.044  0.238  0.647  0.968 

9  0.000  0.016  0.032  0.067  0.310  0.817 

10  0.012  0.016  0.048  0.222  0.635  0.964 

11  0.000  0.015  0.031  0.104  0.340  0.861 

12  0.000  0.012  0.035  0.191  0.572  0.911 

13  0.004  0.019  0.039  0.074  0.377  0.879 

14  0.000  0.008  0.019  0.061  0.290  0.718 

15  0.000  0.016  0.033  0.078  0.393  0.898 

Avg  0.002  0.013  0.033  0.115  0.408  0.846 


1.000  0.7%  0.392  0.131  0.037  0.012  0.000 

1.000  0.832  0.434  0.148  0.047  0.012  0.000 

1.000  0.537  0.183  0.074  0.035  0.008  0.000 

1.000  0.931  0.641  0.252  0.080  0.019  0.000 

1.000  0.540  0.198  0.075  0.036  0.008  0.004 

1.000  0.857  0.480  0.171  0.056  0.012  0.004 

1.000  0.520  0.171  0.063  0.028  0.012  0.000 

1.000  0.807  0.417  0.143  0.046  0.027  0.008 

1.000  0.560  0.198  0.074  0.039  0.012  0.004 

1.000  0.794  0.381  0.125  0.039  0.004  0.000 

1.000  0.939  0.420  0.130  0.057  0.031  0.004 

1.000  0.803  0.377  0.123  0.037  0.012  0.004 

1.000  0.762  0.368  0.129  0.046  0.016  0.003 


Thus,  the  tabulated  normalized  point  response  values 
include  a  component  of  zero-mean  Gaussian  noise  of 
RMS  amplitude: 

0.0013,  midwave 

a  *  (17) 

0.004,  longwave  . 

A  number  of  detector  cells  were  either  totally  inac¬ 


tive  or  too  noisy  to  be  useful,  viz.,  midwave  channels 
14  and  16,  and  longwave  channel  1.  Consequently, 
there  are  neither  point  response  data  nor  Montauk 
Point  clutter  data  corresponding  to  these  channels. 
Although  longwave  channel  16  was  active  at  Montauk 
Point  no  point  response  data  were  available  for  this 
channel.  Finally,  there  were  a  number  of  channels  that 
were  active,  but  for  which  the  point  responses  were 
damaged  by  data  recorder  errors,  viz.,  midwave  chan¬ 
nels  7  and  9,  and  longwave  channel  4. 


A 
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5.  ERROR  PRIMITIVES:  PEAKS  AND  EDGES 


The  input  tc  our  data  correction  algorithm  is  a 
400-sample  time  series, 

x(k)  ,  k  =  1,  2.  ....  400, 

representing  a  single  uncorrected  scan  line  of  1R  data. 
The  elements  of  at(At)  are  integers,  having  the  range 
of  values 

0  2J  x(k)  s  4095. 

The  waveform  x(k)  is  subjected  to  an  initial 
processing  operation  that  generates  as  its  product  a 
vector  of  primitives, 

p(k),  k  =  3,  4,  ....  398, 

having  396  values  in  one-to-one  correspondence  with 
elements  3  to  398  of  x (A)  The  elements  of  p(k)  each 
have  one  of  five  possible  values,  reflecting  an  initial 
judgement  concerning  the  likelihood  of  sample  x(k) 
containing  a  data  error: 

!0  No  initially  obvious  error 
1  Positive  error  peak 
-1  Negative  error  peak  (18) 
2  Positive  error  edge 
-2  Negative  error  edge 

In  this  section  we  explain  the  workings  of  an  al¬ 
gorithm  (PRMTVS)  that  generates  p(k),  implicitly 
defining  the  terms  “error  peak”  and  “error  edge.” 

The  test  for  an  error  peak  at  sample  i  is  conducted 
on  four  sequential  slope  values  (Fig.  20), 

LL  =  x(i  -  1)  -  x(i  —  2) 

L  =  *(/)  -  x(i  -  1) 

R  =  x(i  +  1)  —  x(i)  (19) 

RR  =  x(i  +  2)  -  x(i  +  1). 

The  test  for  an  error  edge  at  sample-pair  (/,  j  +  1) 
is  conducted  on  three  sequential  slope  values  (Fig.  21), 

L  =  x(j)  -  x(j  -  1) 

C  =  x{j  +  1)  -  xV)  (20) 

R  =  x(j  +  2)  -  *0'  +  1)  . 

It  is  in  the  nature  of  our  definitions  for  error  peak 
and  error  edge  that  sequential  elements  of  p(k)  must 
satisfy  the  set  of  adjacency  rules  given  in  Table  5.  For 
example,  it  may  be  inferred  from  Table  S  that  edge 


Figure  20.  Notation  for  error  peak  test.  The  test  for 
an  error  peak  at  sample  I  is  conducted  on  four  sequen¬ 
tial  slope  values  defined  by  the  5-sample  neighborhood 
of  /.  Sample  values  of  x(*)  are  indicated  by  filled  cir¬ 
cles.  Slope  values  are  given  by  LL,  L,  R,  and  RR  (Eq.  19). 


Figure  21.  Notation  for  error  edge  test.  The  test  for 
an  error  edge  at  the  sample-pair  (/,/  +  1)  is  conducted 
on  three  sequential  slope  values  defined  by  the 
4-sample  neighborhood  (/  -  1,  /  +  2).  Sample  values 
of  x{k)  are  indicated  by  filled  circles.  Slope  values  are 
given  by  L.  C,  and  R  (Eq.  20). 


values  must  always  appear  in  groups  of  at  least  two 
sequential  values  (e.g.,  Table  6,  entries  3  and  4),  and 
can  never  appear  in  isolation  (Table  7,  entry  2).  Other 
examples  of  allowed  and  disallowed  primitive  se- 
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quences  are  provided  in  Tables  6  and  7,  respectively. 

The  algorithm  that  extracts  error  peak  and  error 
edge  primitives  from  the  waveform  data  is  called 
PRMTVS  (Appendix  A). 

Tracing  the  logic  flow  in  PRMTVS,  we  see  that  in¬ 
dication  of  an  error  peak  at  sample  k  requires  the  satis¬ 
faction  of  a  necessary  condition,  step  2.11.3;  in 
addition,  one  of  two  auxiliary  conditions  must  also 
be  satisfied,  steps  2.12.1  and  2.13.2. 

Thresholding  the  slope  absolute  value,  step  2.11.3a 
removes  from  consideration  the  very  large  number  of 
low  amplitude  spikes  attributable  to  closed-cover  noise 
(an  additive  Gaussian  noise  having  unity  standard  devi¬ 
ation,  a  =  1).  Step  2.11.3b  is  the  basic  peak  condi¬ 
tion,  assuring  a  sign  change  in  slope. 

Table  5 

Primitive  adjacency  rules.  Table  entry  O  indicates  an 
allowed  adjacency;  entry  X  indicates  a  disallowed  ad¬ 
jacency;  R  indicates  a  required  adjacency. 


p(*)  +1  -1 

+  2  -2 

+  1  X  O 

O  O 

-1  OX 

o  o 

+  2  0  0 

R  O 

-2  0  0 

O  R 

Table  6 

Examples  of  allowed  primitive  sequences. 

P(k) 

Interpretation 

1. 

(...0  1  0...| 

Isolated  peak 

2. 

I...0  1  -1  0...| 

Adjacent  peaks  of 

[...0  1-11  0...) 

alternating  sign 

3. 

[...0  2  2  0...I 

Isolated  edge 

4. 

[...0  2  2  2  2  0...) 

Staircase 

5. 

(...0  2  2  -2  -2  0...I 

Plateau 

6. 

{...0  1  2  2  0...| 

Peak  adjacent  to 
edge 

Table  7 

Examples  of  disallowed  primitive  sequences. 

Adjacency  rule 

p(k) 

(Table  5) 

1.  I...0  I  1  0...|  X 

2.  [...0  2  0... |  R 


The  error  peak  primitives  are  of  two  distinct  types, 
referred  to  by  us  as  PK1  and  PK2,  corresponding  to 
steps  2. 12.1  and  2.13.2,  respectively.  Examples  of  PK1 
and  PK2  error  primitives  are  given  in  Figs.  22  and  23. 

The  quantity  q , ,  defined  by  step  2.13.1  is  inter- 


preted  by  us  as  a  “peak  sharpness  metric”.  Table  8 
provides  values  of  g,  derived  from  the  system  point 
response.  Table  4.  A  threshold  Qt  on  peak  sharp¬ 
ness,  derived  from  Table  8,  is  applied  in  2.13.2  to  as¬ 
sure  that  PRMTVS  does  not  mistakenly  assign  a  value 
p(k)  =  ±  1  to  the  peak  of  an  optical  point  source 


Slope  name  LL  L  R 

Slope  value  -i  4  -4 

Primitive  p<k)  •  0  1 


Figure  22.  Example  of  PK1  error  peak  primitive  at 
sample  /.  The  values  given  for  LL,  L,  and  R  may  be  seen 
to  satisfy  steps  2.11.3  and  2.12.1.  Asterisks  denote  in¬ 
determinate  values  of  p.  In  general,  just  three  slope 
values,  (LL,  L,  and  R,  or  L,  R,  and  RR),  are  needed  to 
determine  a  PK1  error  peak. 


Slope  name  LL  L  R  RR 

Slope  value  2  4  -4-2 

Primitive  pik)  *  0  10 


Figure  23.  Example  of  PK2  error  peak  primitive  at 
sample  i.  The  values  given  for  LL,  L,  R,  and  RR  may 
be  seen  to  satisfy  steps  2.11.3  and  2.13.2.  Asterisks 
denote  indeterminate  values  of  p.  Regarding  2.13.2, 
note  that  for  this  example:  4  =  g,  >  Q,  =  1.85. 
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Table  8 

Metrics  g,  -  <73  are  derived  from  the  system  point  responses.  Values  for  qt  are  ob¬ 
tained  from  step  2.13.1,  Eq.  19,  and  Table  4.  The  value  of  i  to  use  in  Eq.  19  is  i  = 
6  for  midwave  responses  and  /  =  7  for  longwave  responses.  Values  for  q2  are  ob¬ 
tained  from  step  2.25.1,  Eq.  20,  and  Table  4.  Values  for  q3  are  obtained  from  step 
3.32.3  and  Table  4. 


q,  q2  (left)  q2  (right)  <?3 


Channel 

MW 

LW 

MW 

1 

0.593 

* 

0.857 

2 

0.639 

0.750 

0.678 

3 

0.696 

0.726 

0.825 

4 

0.675 

• 

0.826 

5 

0.777 

0.750 

0.866 

6 

0.790 

1.470 

0.769 

7 

• 

1.010 

* 

8 

0.809 

1.448 

0.778 

9 

* 

0.738 

« 

10 

0.894 

1.483 

0.760 

11 

0.921 

0.762 

0.785 

12 

0.936 

1.479 

0.749 

13 

0.995 

0.740 

0.773 

14 

• 

0.778 

• 

15 

1.005 

0.665 

0.717 

16 

• 

• 

* 

minimum 

0.593 

0.665 

0.678 

maximum 

1.005 

1.483 

0.866 

LW 

MW 

LW 

MW 

* 

1.059 

0 

0.578 

0.815 

0.983 

1.506 

0.588 

1.184 

0.983 

0.868 

0.584 

* 

0.953 

* 

0.561 

1.174 

0.941 

0.879 

0.589 

0.735 

0.871 

0.619 

0.586 

1.032 

0 

0.633 

* 

0.730 

0.880 

0.585 

0.579 

1.1% 

• 

0.833 

• 

0.735 

0.832 

0.595 

0.601 

1.392 

0.831 

0.835 

0.608 

0.719 

0.782 

0.641 

0.607 

1.183 

0.734 

0.891 

0.608 

0.836 

• 

1.478 

0 

1.211 

0.697 

0.945 

0.617 

* 

• 

0 

• 

0.719 

0.697 

0.585 

0.561 

1.392 

1.059 

1.506 

0.617 

LW 


0.701 

0.680 
* 

0.692 
0.703 
0.717 
0.685 
0.687 
0.696 
0.709 
0.709 
0.689 
0.713 
0.683 

*  1 

i 

0.680 
0.717 


(such  as  a  small  target  at  long  range)  that  could  pos¬ 
sibly  be  present  in  the  viewed  scene.  Applying  a 
threshold  to  q,  allows  us  to  discriminate  between 
target-like  peaks,  e.g..  Fig.  16,  and  error  peaks,  e.g., 
Fig.  23. 

According  to  Table  8, 

0.59  <  q,  <  1.49, 

from  which  we  surmise  that  setting  a  threshold  value 

q\  <  Q\  =  1-85 

should  assure  that  the  error  peak  criterion  in  PRMTVS 
does  not  inadvertently  cause  damage  to  valid  data,  ex¬ 
cept  on  acceptably  rare  occasions. 

We  momentarily  defer  discussion  of  step  2.13.3  (in¬ 
voking  the  PK2MOD  algorithm). 

Considering  now  the  second  part  of  PRMTVS, 
which  deals  with  error  edges,  we  observe  that  indica¬ 


tion  of  an  error  edge  at  sample  k  requires  the  satis¬ 
faction  of  a  necessary  condition,  step  2.21.2.  In 
addition,  auxiliary  conditions  must  also  be  satisfied: 
either  2.22.1,  2.23.1,  or  2.25.2. 

Thresholding  the  slope  absolute  amplitude,  2.21 .2a, 
removes  most  noise  artifacts  from  consideration.  Step 
2. 2 1. 2b  imposes  the  requirement  that  neither  candi¬ 
date  edge  point  has  previously  been  labeled  by 
PRMTVS  as  an  error  peak. 

The  error  edge  primitives  are  of  two  distinct  types, 
called  EDG1  and  EDG2,  corresponding  respectively 
to  steps  2.22.1/2.23.1,  and  2.25.2.  Examples  of  EDG1 
and  EDG2  error  primitives  are  given  in  Figs.  24  and  25. 

The  quantity  q2  defined  in  step  2.25.1  is  interpret¬ 
ed  by  us  as  an  “edge  sharpness  metric.”  Table  8  pro¬ 
vides  values  of  q2  derived  from  the  system  point 
responses.  Table  4.  A  threshold,  Q2 ,  on  edge  sharp¬ 
ness,  derived  from  Table  8,  is  applied  in  step  2.25.2 
to  assure  that  PRMTVS  does  not  mistakenly  assign 
an  error  value p(  k)  =  ±  2  to  the  shoulder  of  an  opti¬ 
cal  point  source  possibly  present  in  the  viewed  scene. 
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Sample,*  /-I 

Slope  name 
Slope  value 
Primitive  ptk)  * 


Figure  24.  Example  of  EDG1  error  edge  primitive  at 
sample  pair  (/,  i  +  1).  The  values  given  for  L,  C,  and 
R  may  be  seen  to  satisfy  steps  2.21.2  and  2.22.1.  Aster¬ 
isks  denote  indeterminate  values  of  p. 


Sample,  *  <  -  1  /  /  +  1  /  +  2 

Slope  name  L  C  R 

Slope  value  1  6  1 

Primitive  ptk)  *  2  2 

Figure  25.  Example  of  EDG2  error  edge  primitive  at 
sample  pair  (/,  /  +  1).  The  values  given  for  L,  C,  and 
R  may  be  seen  to  satisfy  steps  2.21.2  and  2.25.2.  Aster¬ 
isks  denote  indeterminate  values  of  p.  Regarding 
2.25.2,  note  that  for  this  example:  3  =  q2  >  Q2  =2. 


The  error  primitives  defined  above  are  unable  to  de¬ 
tect  a  type  of  length-three  artifact  illustrated  in  Fig. 
26.  (The  likelihood  of  this  artifact  type  assembling  by 
chance,  i.e.,  frequency-of-occurrence  statistics,  has  not 
been  established.)  The  PRMTVS  algorithm  at  step 
2.13.3  invokes  a  routine,  PK2MOD  (cf.  Appendix  A), 
that  augments  our  original  definition  of  the  PK2  er¬ 
ror  primitive  by  implementing  a  six-slope  test  capa¬ 
ble  of  discovering  “flat-footed  peak”  artifacts  such 
as  that  depicted  in  Fig.  26.  However,  we  note  that 
PK2MOD  was  not  operational  during  the  numerical 
validation  studies  discussed  in  connection  with  Figs. 
12  through  14. 


Sample  number:  *  1  2  3  4  5  6  7 

Original  data:  a(k)  0  1  2  3  2  1  0 

Errors:  s(k)  0  0  9  14  9  0  0 

Damaged  data:  b(k)  0  1  11  17  11  1  0 

PK2:  p(k)  0  0  0  0  0  0  0 

PK2MOD:  p(k)  0  0  0  10  0  0 

Figure  26.  Example  of  a  length-3  artifact  that  would 
go  undetected  but  for  the  action  of  algorithm  PK2MOD. 
Errors  of  amplitude  s(k)  are  implanted  at  *  =  3, 4,  and 
5.  For  k  =  5,  block  2.25  calculates  g2  =  1.43  <  02 
=  2;  i.e.,  the  shoulders  of  the  artifact  are  not  steep 
enough  to  qualify  as  error  edges.  Fork  =  4,  block  2.13 
calculates  g,  =  1.2  <  Q,  =  1.85;  i.e.,  the  4-slope  PK2 
test  also  fails  to  discover  the  artifact.  However,  the 
6-slope  test  implemented  in  PK2MOD  returns  a  non¬ 
null  value  for  p(4),  as  desired.  As  determined  by 
PK2MOD,  the  tails  of  the  peak  depicted  here  flatten 
too  quickly  (compare  with  Fig.  15a). 

At  step  3.32.4  in  PK2MOD  the  shape-related  met¬ 
ric  q,  is  compared  with  a  threshold  value  Q, .  The 
value  assigned  to  Qy  is  developed  by  applying  the  six- 
slope  test  in  PK2MOD  to  point  response  data  on  a 
channel-by-channel  basis  (cf.  Table  8),  in  a  process 
strictly  analogous  to  the  development  of  threshold 
values  for  q,  and  q:,  the  peak  and  edge  sharpness 
metrics  defined  earlier. 
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6.  INITIAL  ESTIMATES  FOR  ARTIFACT  BOUNDARIES 


The  vector  of  primitives  generated  by  the  PRMTVS 
algorithm,  p(k),  is  next  operated  on  by  the  BOUNDS 
algorithm  to  produce  estimates  for  the  starting  and 
ending  locations  of  each  error  artifact  in  the  scan  line. 
We  define  the  beginning  sample  number  of  artifact 
j,  b(j),  as  the  first  defective  datum  in  artifact  j  (not 
the  last  valid  datum  preceding  the  artifact).  Similar¬ 
ly,  the  final  sample  number  of  artifact  j,  f(j ),  we  de¬ 
fine  as  the  last  defective  datum  in  artifact  j  (not  the 
first  valid  datum  following  the  artifact). 

The  BOUNDS  algorithm  (Appendix  A)  begins  by 
searching  for  the  first  non-null  element  oip(k)  (step 
4.11.2).  If  a  non-null  element  is  found  (say,  at  loca¬ 
tion  kt ),  the  artifact  counter  is  incremented  to  unity 
(n  =  1)  and  the  beginning  location  of  the  first  artifact, 
/HD,  is  estimated  (step  4.12.2)  as  /HD  =  k, .  The  al¬ 
gorithm  then  looks  for  three  successive  null  values  of 
p(k)  before  deciding  that  the  artifact  has  terminated 
(step  4.13.2). 


The  operation  of  PRMTVS  and  BOUNDS  is  illus¬ 
trated  in  Fig.  27  for  an  isolated  error  implanted  into 
a  constant-valued  waveform.  For  this  simple  exam¬ 
ple  the  artifact  boundary  estimates  generated  by 
BOUNDS  precisely  indicate  the  error’s  location. 

We  define  a  “plateau”  as  an  error  artifact  having 
initial  and  final  primitives  that  are  edges  of  opposite 
sign,  i.e.,  either  (p(/>),p(/»  =  (2,  -  2),  or  [p(b),p(J)\ 
=  (-2,2),  so  that  p(b)  x  p(f)  =  -4  (step  4. 12.8). 
Our  estimates  for  the  boundaries  of  a  plateau  artifact 
involve  a  shrinking  process  (block  4.14)  rationalized 
in  Fig.  28. 

We  define  a  first  estimate  for  the  length  of  artifact 
j,  i.e.,  the  number  of  sequential  defective  samples  in 
the  artifact  as 


L(J)  =  fU)  ~  b(j)  +  1  • 


Sample  no. :  k  1 
Errors'.  stk)  0 
PRMTVS:  p(k)  0 

BOUNDS: 


5  6  7 
0  0  0 
0  0  0 


Figure  27.  Isolated  error  implanted  at  location  k  =  4  in  a  constant-valued  waveform. 
The  PRMTVS  algorithm  generates  a  peak  error  primitive  at  the  location  of  the  error,  i.e., 
p(4)  =  1,  and  null  primitives  elsewhere.  The  BOUNDS  algorithm  indicates  the  existance 
of  an  error  artifact  having  beginning  and  final  locations  (b  and  /,  respectively)  coincident 
with  one  another  and  with  the  error’s  location. 


Sample  no.:  k  3 
Errors:  stk)  0 

PRMTVS:  ptk)  0 

BOUNDS: 


8  9  10  11 

0  0  0  0 

-2  0  0  0 


Figure  28.  Plateau  error  artifact  implanted  at  locations  k  =  6  and  7  in  a  constant-valued 
waveform.  The  BOUNDS  algorithm  executes  a  boundary-shrinking  operation  for  this  type 
of  artifact  such  that  the  estimated  beginning  and  final  artifact  locations  (b  and  f,  respec¬ 
tively)  coincide  with  the  locations  of  the  implanted  errors  (block  4.14).  A  plateau  is  de¬ 
fined  as  an  artifact  for  which  lp(b)’P(/)]  =  -4  (step  4.12.8). 
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7.  IDENTIFICATION  OF  ISOLATED  ERRORS 


Contrary  to  the  impression  perhaps  conveyed  by 
Figs.  27  and  28,  the  artifact  boundary  estimates  gener¬ 
ated  by  BOUNDS  often  provide  an  inaccurate  indi¬ 
cation  of  precisely  which  data  are  in  error.  This  point 
is  illustrated  by  analysis  of  five  examples  (Figs.  29 
through  33).  The  lessons  learned  from  these  examples 
are  generalized  to  a  set  of  rules  for  refining  the  initial 
artifact  boundary  estimates. 

Figure  29  depicts  a  single  error  implanted  into  a 
constant  slope  waveform.  The  slope,  though  greater 
than  3  (step  2. 1 1.3a),  is  much  smaller  than  the  ampli¬ 
tude  of  the  implanted  error.  The  simple  peak-finding 
algorithm  in  PRMTVS  has  labeled  sample  k  =  5  in 
Fig.  29  as  a  negative  peak,  and  BOUNDS  is  designed 
to  retain  all  such  labeled  features  as  candidate  errors, 
i.e.,  as  candidates  for  correction.  In  this  instance,  the 
artifact  boundary  estimates  provided  by  BOUNDS  en¬ 
compass  the  actual  error  at  k  =  4  and  a  valid  datum 
at  k=  5. 


Sample  no. :  *  1234  5678 

Error:  s(k)  0  0  0  e0  0  000 

PRMTVS:  pfk)  0  1-1  0  0  0 

BOUNDS:  b  f 

Figure  29.  A  single  error  is  implanted  at  location  k 
=  4  in  a  constant-slope  waveform.  The  slope  ampli¬ 
tude  is  greater  than  3,  but  much  smaller  than  e0  (the 
amplitude  of  the  implanted  error).  The  artifact  bound¬ 
ary  estimates  provided  by  BOUNDS  encompass  the  er¬ 
ror  at  k  =  4  and  a  valid  datum  at  k  =  5. 


Figure  30  is  similar  to  Fig.  29  in  once  again  depict¬ 
ing  a  single  error  implanted  into  a  constant-slope  wave¬ 
form.  However,  in  Fig.  30  the  amplitude  of  the  im¬ 
planted  error  equals,  by  coincidence,  the  slope  am¬ 
plitude.  Once  again,  the  boundaries  provided  by 
BOUNDS  encompass  the  actual  error  at  k  =  4. 
However,  three  valid  data  points  (k  =  3,  5,  and  6) 
are  also  indicated  as  candidates  for  correction,  indicat¬ 


ing  that  BOUNDS  is  apparently  unable  to  localize  the 
error. 

Figure  31  depicts  the  implantation  of  a  single  error 
at  the  peak  of  the  point  response,  giving  rise  to  a 


Sample  no.:  *123456789 
Error:  s(k)  000  e0  00000 

PRMTVS:  p(k)  2  2  2  2  0  0  0 

BOUNDS:  b  f 

Figure  30.  A  single  error  is  implanted  at  location  k 
=  4  in  a  constant-slope  waveform.  The  amplitude  of 
the  implanted  error  equals  by  coincidence  the  slope 
amplitude.  The  artifact  boundary  estimates  generated 
by  BOUNDS  encompass  the  error  at  k  =  4  and  three 
valid  data  points  as  well. 


Sample  no. :  *123456789 
Error:  sfk)  0  0  0  0  -Sq  0  0  0  0 

PRMTVS:  pfk )  0  1-1  10  0  0 

BOUNDS:  b  f 

Figure  31.  A  single  error  occurs  at  the  peak  of  the 
point  response,  giving  rise  to  a  double-peaked  artifact. 
Artifact  boundaries  b  and  /  estimated  by  BOUNDS  en¬ 
compass  the  error  at  k  =  5,  as  well  as  valid  data  at 
k  =  4  and  k  =  6. 
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Sample  no.:  A  1  2  3 

Error:  s(k)  0  0  0 

PRMTVS:  p(k)  0 

BOUNDS: 


0  0 
2  2  2 


8  9  10 
0  0  0 
0  0  0 


Figure  32.  A  single  error  occurs  one  sample  distant 
from  the  peak  of  a  flat-topped  point  response  (Fig.  15b), 
creating  a  staircase-shaped  artifact.  BOUNDS  labels 
a  4-sample  sequence  as  suspect. 


Sample  no.:  *  123  4  56789 

Error:  s(k)  000  -e0  00000 

PRMTVS.  p/kf  0  0  1  0  0  0  0 

BOUNDS:  b/f 

Figure  33.  A  single  error  occurs  one  sample  distant 
from  the  peak  of  an  isosceles  point  response  (Fig.  15a). 
The  boundary  estimates  b  and  f  provided  by  BOUNDS 
do  not  encompass  the  error  at  k  =  4  for  the  depicted 
example. 


double-peaked  artifact.  The  artifact  boundaries  gener¬ 
ated  by  BOUNDS  encompass  the  error  at  k  =  5,  as 
well  as  one  valid  datum  on  each  side  of  the  error. 

In  Fig.  32  a  single  error  occurs  one  sample  distant 
from  the  peak  of  a  flat-topped  point  response  (e.g., 
Figs.  13b  and  17).  The  amplitude  of  the  error  is  such 


that  a  staircase  artifact  (Table  6,  entry  4)  is  created. 
The  four-sample  sequence  labeled  as  suspect  by 
BOUNDS  encompasses  the  isolated  error. 

As  a  final  canonical  example  of  single  sample 
damage  we  consider  Fig.  33,  in  which  the  single  error 
occurs  one  sample  distant  from  the  peak  of  an  isosceles 
point  response  (e.g..  Figs.  15a  and  16).  The  scale  of 
the  point  response  and  the  amplitude,  sign,  and  loca¬ 
tion  of  the  error  all  conspire  to  produce  coincident 
boundary  estimates  b  and  /that,  in  this  pathological 
instance,  do  not  encompass  the  actual  error. 

We  draw  two  conclusions  from  the  preceding  ex¬ 
amples,  relevant  to  the  design  of  the  ISOERR  al¬ 
gorithm,  whose  purpose  is  to  detect  and  correct 
isolated  errors. 

•  An  artifact  having  a  BOUNDS  length  estimate 
as  large  as  4  may  be  caused  by  an  isolated  error 
(e.g..  Figs.  30  and  32).  Thus,  a  single  sample  cor¬ 
rection  should  be  attempted  whenever  L  <  4. 

•  An  artifact  having  a  BOUNDS  length  estimate 
of  1  may  be  caused  by  an  isolated  error  located 
one  sample  distant  (left  or  right)  from  the  non¬ 
null  primitive  (Fig.  33). 

The  ISOERR  algorithm  (Appendix  A)  begins  with 
an  initial  determination  of  whether  a  single  sample  cor¬ 
rection  should  be  attempted  for  artifact  j.  As  noted 
above,  single  sample  correction  is  attempted  whenever 
L(j)  s  4  (step  5.11.2).  If  L(j)  >  4,  the  artifact  is 
left  untreated  by  ISOERR. 

Quantities  B  and  F,  initialized  in  step  5.12.1,  are 
the  boundaries  of  a  search  interval:  each  sample  value 
of  x(n)  (for  B  <  n  <  F)  is  separately  considered  as 
a  possible  isolated  error. 

As  discussed  in  connection  with  Fig.  33,  when  L  ( j ) 
=  1  we  must  consider  the  possibility  that  an  isolated 
error  may  have  occurred  at  any  of  three  possible  lo¬ 
cations,  viz.,  b  -  1,  6,  and  b  +  1.  Thus,  ISOERR 
next  tests  to  see  if  L  =  1  (step  5.12.2);  if  so,  ISOERR 
expands  in  block  5.125  the  search  boundaries  B  and  F. 

We  now  draw  what  we  refer  to  as  the  “block  5.13 
hypothesis.”  In  block  5. 13  we  test  the  hypothesis  that 
the  yth  artifact  is  caused  by  an  isolated  error  occur¬ 
ring  at  the  location  n  within  the  artifact. 

The  first  action  within  block  5. 1 3  is  to  save  the  origi¬ 
nal,  presumably  defective,  data  value  at  location  n 
(step  5.13.1).  Next ,  the  NTRPL8  algorithm  is  used  in 
step  5.13.2  to  generate  a  replacement  value  for  x(n), 
called  y(n),  by  interpolation  from  the  neighboring 
values  x(n  -  2),  x{n  -  1 ),  x(n  +  1 ),  and  x(n  + 
2),  which  are  hypothesized  to  be  valid  data.  (The  cu¬ 
bic  interpolation  implemented  by  NTRPL8  is  derived 
in  Appendix  B).  If  the  block  5.13  hypothesis  is  via¬ 
ble,  the  artifact  will  now  have  been  eliminated.  Thus 
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where  we  recall  in  applying  Eq.  22  that  waveform  x(k ) 
has  been  modified  by  replacing  jt[A/(i)J  with  the  in¬ 
terpolated  value  y[JV(/)j  (step  5.13.3). 

For  example,  considering  Fig.  34,  the  first  of  the 
viable  single  sample  corrections  (i  -  1 )  occurs  at  lo¬ 
cation  k  =  N(\)  =  3;  the  second  viable  single  sample 
correction  (/  =  2)  occurs  at  location  k  =  N(2)  =  4. 

We  close  this  section  with  one  more  example  illus¬ 
trating  the  interaction  between  the  algorithms  dis¬ 
cussed  up  to  this  point:  PRMTVS,  BOUNDS, 
ISOERR,  and  NTRPL8. 

Figure  35  depicts  a  constant-valued  waveform  into 
which  equal  amplitude  errors  have  been  inserted  at  lo¬ 
cations  k  =  4  and  k  =  6.  The  damaged  waveform 
is  input  to  PRMTVS,  which  generates  a  vector  of 
primitives, /?(£),  containing  three  non-null  characters 
(Fig.  35).  The  BOUNDS  algorithm  identifies  one  er¬ 
ror  artifact  in  the  waveform,  extending  from  b  =  4 
to/  =  6,  encompassing  the  data  errors  at  k  =  4  and 
k  =  6,  and  a  valid  datum  at  k  =  5.  Since  L  =  3  s 
4,  ISOERR  attempts  to  eliminate  the  artifact  with  a 
single  sample  correction  applied  at  each  of  three  lo¬ 
cations:  k  -  4,  5,  and  6. 


Sample  no.:  k 
Errors:  stk) 

PRMTVS:  p/k) 


»0 

1 


0 

-1 


e0 

1 


Tarwsrmrw 


On  the  first  pass  through  block  5.13  (w  =  4), 
ISOERR  applies  a  trial  correction  at  location  k  =  4. 


As  indicated  in  Fig.  35,  subsequent  output  generated 
by  PRMTVS  at  step  5.13.4,  p'(Ar),  has  a  non-null 
character,  p' (6)  =  1 .  Consequently,  block  5.135  is  not 
executed,  and  the  counter  of  viable  corrections  retains 
its  initial  value,  t  =  0. 


On  the  second  pass  through  block  5.13  ( n  =  5), 
ISOERR  applies  a  trial  single  correction  at  location 
k  =  5.  The  replacement  value  for  x(5)  generated  by 
NTRPL8  is 


y( 5)  =  *(5)  +  4e0/3 


as  indicated  in  Fig.  35.  Inspection  of  Fig.  35  shows 
that  the  originally  separate  errors  at  k  =  4  and  k  = 
6  have  been  “blended”  smoothly  with  yet  a  third  er¬ 
ror  produced  by  ISOERR  itself.  Indeed,  when 
PRMTVS  is  applied  to  the  resultant  waveform  (now 
containing  three  errors),  it  is  found  that  the  negative 
error  peak  originally  present  at  k  =  5  has  been  elimi¬ 
nated,  i.e.,  while  originally  p( 5)  =  -  1,  after  the  se¬ 
cond  pass  through  block  5.13  we  obtain  p'(5)  -  0. 
However,  the  “shoulders”  of  the  error  artifact  are  still 
too  abrupt  to  be  optically  induced:  two  error  edges 
now  appear  in  p'  (k). 


It  follows  from  the  foregoing  discussion  that 
ISOERR  leaves  unchanged  the  multiple  error  artifact 
depicted  in  Fig.  35:  ISOERR  neither  damages  good 
data  n or  corrects  either  of  the  errors. 


i 

3 


3 


Block  5.13 

First  pass:  p'lk)  0  0  0  1  0  0 

Second  pass:  p'lk)  2  2  0  -2  -2  0 

Third  pass:  p’lk)  0  1  0  0  0  0 

Figure  35.  Equal  amplitude  errors  are  inserted  at  lo¬ 
cations  k  =  4  and  k  =  6  into  a  constant-valued  wave¬ 
form.  The  dashed  waveform  depicted  here,  passing 
through  y( 5),  is  the  waveform’s  appearance  after  the 
second  pass  through  step  5.13.3. 


We  note  that  ISOERR  detects  and  corrects  with  high 
accuracy  all  six  examples  of  single-sample  damage 
presented  above  in  Figs.  29  through  34. 


Finally,  we  note  that  the  version  of  EDS  used  in  our 
statistical  validation  study  had  a  call  to  PK3MOD  fol¬ 
lowing  step  5.13.4,  rather  than  at  step  2.13.3,  where 
PK3MOD  could  potentially  have  had  greater  impact 
in  reducing  undetected  errors. 
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8.  IDENTIFICATION  OF  MULTIPLE  ERRORS 


Multiple  sequential  errors  generally  cannot  be  cor¬ 
rected,  i.e.,  there  exists  no  way  to  reconstruct  with  high 
accuracy  the  original  data  values.  The  most  we  can 
accomplish  in  our  treatment  of  multiple  errors  is  to 
identify  which  of  the  data  are  defective,  and  then  to 
smoothly  interpolate  through  the  damaged  regions. 
Once  data  errors  are  properly  labeled  they  may  be  ex¬ 
cluded  from  use  in  subsequent  analysis  (“excised” 
from  the  data  set).  Consequently,  the  net  effect  of  hav¬ 
ing  multiple  errors  is  a  speckling  of  “blind  spots”  in 
the  data,  with  a  consequent  reduction  in  data  set  size 
and  ease  of  use,  but  no  degradation  in  the  structural 
integrity  of  the  waveforms. 

Our  algorithm  for  detecting  and  smoothing  multi¬ 
ple  sequential  errors  is  called  MLTPLS  (Appendix  A). 
The  waveform  data  entered  as  input  to  MLTPLS  have 
previously  been  operated  on  by  ISOERR;  i.e.,  the  iso¬ 
lated  errors  have  been  corrected  prior  to  operating  on 
the  data  with  MLTPLS. 

In  step  7.12.2  MLTPLS  attempts  to  excise  artifact 
n  by  linearly  interpolating  from  [B(n)  -  1]  to  (F(n) 
+  1].  The  smoothed  version  of  x{k)  is  then  examined 
by  PRMTVS  in  step  7.12.3  to  determine  whether  any 
excessively  abrupt  peaks  or  edges  remain  in  the  wave¬ 
form.  If  so,  the  artifact  boundaries  are  expanded  in 
block  7.31  to  obtain  new  values  for  B(n)  and  F(n), 
and  MLTPLS  iterates  once  again  through  the  linear 
interpolation  block.  The  process  of  boundary  expan¬ 
sion  and  artifact  smoothing  continues  until  PRMTVS 
returns  a  vector  of  null  primitives,  which  causes 
MLTPLS  to  terminate. 

It  was  noted  during  algorithm  validation  that  many 
of  the  undetected  errors  (Fig.  13)  occur  in  between 
pairs  of  recognized  multiple  error  artifacts.  This  ob¬ 
servation  motivated  the  development  of  the  MERG¬ 
ER  algorithm  (Appendix  A),  invoked  from  MLTPLS 
at  step  7.31.3.  However,  our  statistical  validation, 
(Figs.  12  through  14)  is  based  on  a  version  of  EDS 
that  does  not  include  MERGER. 

The  values  assigned  to  constants  C,  and  C2  at  step 
9.11.1  establish  a  trade-off  between  data  damage  and 
missed  errors.  Quantifying  this  trade-off,  and  thus 
selecting  optimum  values  for  C,  and  C2 ,  would  re¬ 
quire  an  extended  Monte  Carlo  analysis  that  has  not 
been  performed. 

We  next  illustrate  the  operation  of  MLTPLS  with 
three  examples. 

Figure  36  depicts  a  large  error  artifact  implanted 
in  a  slowly  varying  waveform.  In  this  instance,  the 
original  boundary  estimates  b  and  /  generated  by 


Sample  no.: 
Errors: 
PRMTVS: 
BOUNDS: 
Step  7.12.3 


/  _ _  .  ■  . . 

k  3  4  5  6  7  8  9  10  11  12  1314  15  16 

s(k)  0  0  0  e,  0  83  e4~e5  eg  0  0  0  0 

p(k!  0  0  0  1  -1  0  0  1  -1  1  -1  0  0  0 

b  ( 

p'(k)  0000  000  0000000 


Figure  36.  Example  of  a  large  error  artifact  requiring 
for  excision  just  one  pass  through  block  7.12.  As 
shown  here,  linear  interpolation  through  the  excised 
region  produces  small  correction  errors  when  the  origi¬ 
nal  data  (open  circles)  are  slowly  varying.  Step  1.7  in 
executive  routine  EDS  restores  the  measured  value 
when  it  differs  only  slightly  from  the  corrected  value. 
This  restoration  step  prevents  damage  to  samples  x(8) 
and  x(13)  in  the  depicted  example. 


BOUNDS  do  not  need  to  be  modified,  and  only  one 
pass  through  block  7.12  is  required.  Moreover,  the 
differences  between  the  actual  data  values  (O)  and  the 
interpolated  values  (X)  (i.e.,  the  “correction  errors”) 
are  apparently  small. 

In  Fig.  37  two  sequential  data  errors  are  implanted 
in  an  otherwise  constant-valued  waveform.  In  this  in¬ 
stance,  the  very  large  error  at  k  =  7  prevents  PRMTVS 
from  recognizing  the  smaller  error  at  k  =  8  either  as 
an  error  peak  or  an  error  edge.  However,  after  one 
pass  through  block  7.12  has  largely  suppressed  the 
dominant  error  at  k  =  7,  PRMTVS  on  the  second  pass 
through  block  7.12  successfully  uncovers  the  smaller 
error  at  k  =  8. 

in  Fig.  38  we  illustrate  the  point  that  multiple  er¬ 
rors  occurring  in  a  high-curvature  waveform  typical¬ 
ly  produce  significant  irrecoverable  damage.  In  this 
figure,  errors  occur  at  two  locations  in  the  center  of 
the  longwave/channel  2  point  response.  On  the  first 
pass  through  block  7.12,  linear  interpolation  from  ( B 
-  I)t  to(F+  1),  produces  the  dashed  curve  shown 
in  Fig.  38.  On  successive  passes,  the  residual  highly 
damaged  structure  is  progressively  sliced  away  until, 
on  the  third  pass,  MLTPLS  terminates  with  the  data 
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Sample  no.: 
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10 

11 

Errors: 

s(k) 

0 

0 

0 

3e0 

0 

0 

0 

PRMTVS. 

plk) 

0 

0 

0 

1 

0 

0 

0 

0 

BOUNDS. 

b/f 

Step  7.12.3 
First  pass: 

p'(k) 

0 

0 

2 

2 

1 

0 

0 

0 

Second  pass: 

P'lkt 

0 

0 

0 

0 

0 

0 

0 

0 

Sample  no.: 
Errors: 
PRMTVS: 
BOUNDS: 


<F+1), 


p(k)  0  0 


1  -1 
b 


Figure  37.  Example  of  a  length-2  artifact  requiring  for 
excision  two  passes  through  block  7.12.  Errors  are  im¬ 
planted  at  two  locations,  k  =  7  and  k  =  8,  in  a  con¬ 
stant  amplitude  waveform.  On  the  first  pass,  the 
SMOOTH  algorithm  performs  linear  interpolation  from 
(8  -  1)1  to  (F  +  1), .  As  shown  by  the  dashed  wave¬ 
form,  the  artifact  has  not  yet  been  eliminated.  On  the 
second  pass  through  7.12,  SMOOTH  performs  linear 
interpolation  from  (8  -  1)2  to  (F  +  i)2.  Algorithm 
MLTPLS  then  terminates  with  the  artifact  having  been 
wholly  excised. 


Step  7.12.3 

First  pass:  p'(k)  00  000  0100 

Second  pass:  p'(k)  02  200  0000 

Third  pass:  p'(k)  00  000  0000 

Figure  38.  Errors  occur  at  two  locations  in  the  center 
of  a  point  response  waveform.  On  the  first  pass 
through  block  7.12,  linear  interpolation  from  (8  -  1), 
to  (F  +  1),  produces  the  dashed  curve  shown  in  the 
figure.  Algorithm  MLTPLS  terminates  after  the  third  in¬ 
terpolation,  from  sample  (8  -  l)3  to  sample  (F  +  l)3. 
Thus,  a  5-sample  blind  spot  is  created,  masking  out 
both  the  data  errors  and  the  point  response. 


errors  and  point  response  both  having  been  wholly  ex¬ 
cised.  The  final  interpolation  is  between  samples  (B 
—  1  )3  =6  and  ( F  +  1)3  =  12.  Thus,  a  5-sample 
blind  spot  has  been  created  in  the  data,  masking  out 


what  had  originally  been  a  significant  structural 
feature. 

The  peak-slicing  process  illustrated  in  Fig.  38  is  also 
illustrated  in  Figs.  6,  8,  and  1 1  (feature  D). 


9.  EXECUTIVE  ROUTINE  EDS 

i 


Now  that  the  elements  of  our  error  detection  sys¬ 
tem  have  been  discussed,  we  bring  all  of  our  algorithms 
together  under  the  control  of  executive  routine  EDS 
(Appendix  A). 

In  step  1.1,  EDS  uses  PRMTVS  to  screen  a  scan 
line  of  data.  If  PRMTVS  detects  no  excessively  abrupt 
peaks  or  edges  in  the  data,  EDS  terminates  in  step  1 .2. 
Otherwise,  the  locations  of  the  detected  error  peaks 


and  error  edges  are  used  by  BOUNDS  in  step  1.3  to 
determine  neighborhoods  within  which  error  smooth¬ 
ing  will  be  required. 

In  step  1 .4,  EDS  applies,  via  1SOERR,  a  single  sam¬ 
ple  correction  to  each  of  the  damaged  neighborhoods 
amenable  to  this  simplest  expedient.  Assuming  that 
at  least  one  neighborhood  has  been  repaired  with  a 
single  sample  fix,  block  1.51  operates  on  the  partially 
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repaired  waveform  to  establish  a  residual  set  of  mul¬ 
tiply  damaged  neighborhoods.  The  multiple  error 
regions  are  smoothed  by  the  MLTPLS  algorithm  in 
step  1.6. 

The  last  operation  in  EDS,  step  1.7,  compares  the 
smoothed  waveform  with  a  stored  version  of  the  origi¬ 
nal  unprocessed  waveform.  Where  data  values  have 


been  changed  by  only  a  small  amount  (nominally,  four 
or  fewer  digital  counts)  the  original  data  value  is  re¬ 
stored  to  the  waveform.  (Sensor  internal  noise  is  one 
digital  count,  RMS;  sensor  full-scale  dynamic  range 
is  40%  digital  counts.) 

We  note  that  step  1 .7  was  inactive  during  the  statisti¬ 
cal  validation  discussed  in  connection  with  Figs.  12 
through  14. 


10.  CONCLUSIONS 


A  Navy  program  was  recently  begun  to  obtain  a  test 
set  of  infrared  (IR)  background  images  for  develop¬ 
ing  IR  search  and  track  clutter  suppression  techniques. 
Unfortunately,  the  initial  set  of  measurements  was 
contaminated  by  a  malfunctioning  digital  recorder 
with  about  4%  defective  data  (3%  small-amplitude  er¬ 
rors  and  \°Io  large-amplitude  errors). 

In  this  report  we  present  a  new  Error  Detection  and 
Smoothing  algorithm,  EDS,  and  we  demonstrate  that 
EDS  is  capable  of  eliminating  from  the  IR  data  almost 
all  large  amplitude  errors  (amplitude  greater  than  9 
times  the  standard  deviation  of  the  sensor  closed-cover 
noise). 

The  EDS  algorithm  is  generally  capable  of  dis¬ 
criminating  electrically  induced  burst  noise  (including 
burst  noise  from  electromagnetic  interference)  from 
structural  features  induced  by  optical  variations  in  the 
viewed  scene.  Moreover,  EDS  is  generally  applicable 
to  data  from  any  optical  scanner  whose  output  is  sam¬ 
pled  digitally  multiple  times  per  dwell. 

The  EDS  algorithm  operates  entirely  in  the  time  do¬ 
main.  The  bases  for  error  detection  are  smoothness 
criteria  extracted  from  measurements  made  against  a 
laboratory  point  source,  and  do  not  depend  on  the 
highly  variable  structure  of  the  IR  clutter  data  viewed 
in  the  field.  In  effect,  EDS  exploits  the  redundancy 
inherent  in  highly  sampled  image  data. 

The  computation  is  structured  such  that  the  data 
are  pre-screened  with  a  simple  algorithm  that  allows 
further  attention  to  be  concentrated  only  on  “suspect” 
neighborhoods.  There  are  no  computation-intensive 
steps,  such  as  FFT’s  or  histograms,  required  for  the 
data  set  as  a  whole. 

The  data  are  sampled  at  a  sufficiently  high  rate  that 
isolated  errors  can  be  corrected  by  EDS  to  high  ac¬ 


curacy.  However,  multiple  sequential  errors  general¬ 
ly  cannot  be  corrected;  i.e.,  there  exists  no  way  to 
reconstruct  with  high  accuracy  the  original  data  values. 
Our  approach  to  treating  multiple  errors  is  to  identi¬ 
fy  which  of  the  data  are  defective,  and  then  to  smooth¬ 
ly  interpolate  through  the  damaged  regions.  Once  data 
errors  are  properly  labeled  they  may  be  excluded  from 
use  in  subsequent  analysis  (“excised”  from  the  data 
set).  Consequently,  the  net  effect  of  having  multiple 
errors  is  a  speckling  of  “blind  spots”  in  the  data,  with 
a  consequent  reduction  in  data  set  size  and  ease  of  use, 
but  no  degradation  in  the  structural  integrity  of  the 
waveforms. 

A  valuable  feature  of  the  new  algorithm  is  its  abili¬ 
ty  to  accurately  isolate  error  bursts  (i.e.,  long  sequences 
of  defective  data)  by  means  of  a  region  growing  ap¬ 
proach. 

We  have  presented  a  statistical  validation  based  on 
evaluating  the  performance  of  EDS  against  a  synthet¬ 
ically  degraded  data  set.  We  have  found  that,  on  the 
average,  less  than  1%  of  the  iarger-amplitude  errors 
propagate  undetected  through  the  error  detection  al¬ 
gorithm.  The  density  of  large  amplitude  errors  is  thus 
reduced  from  about  60  per  frame  to  about  0.5  per 
frame  (counting  only  those  errors  whose  amplitude  is 
at  least  10  times  the  sensor  noise). 

We  have  observed  that  the  errors  undetected  by  EDS 
tend  to  be  of  the  following  kinds: 

•  Errors  at  the  edge  of  a  recognized  error  burst, 
where  the  error  burst  abuts  a  region  of  high  cur¬ 
vature  in  the  IR  data, 

•  Errors  sandwiched  closely  between  a  pair  of 
recognized  error  artifacts  (typically,  of  large 
extent), 
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•  Errors  occurring  in  waveforms  that  are  identified 
by  EDS  as  having  exceptionally  high  error 
content. 

These  observations  could  be  used  to  derive  useage 
guidelines  for  the  data.  Formalized  as  algorithms  (e.g., 
the  MERGER  algorithm  in  Appendix  A),  it  is  proba¬ 
ble  that  these  guidelines  could  reduce  the  rate  of  un¬ 
detected  errors  significantly  below  the  nominal  0.5 
errors  per  frame  quoted  above. 

Finally,  it  must  be  said  that  the  block-structured  al¬ 
gorithms  presented  as  EDS  in  Appendix  A  are  com¬ 
pletely  restructured  from  the  FORTRAN  source  that 


was  used  in  our  statistical  validation.  It  is  possible  that, 
in  the  process  of  casting  EDS  into  a  more  easily  ex¬ 
plained  form,  we  may  have  introduced  logic  errors. 

We  conclude  that  our  immediate  objectives  were  met 
in  the  sense  that  EDS  proved  adequate  for  extracting 
from  the  damaged  IR  data  a  limited  interim  data  set 
for  use  until  better  data  were  obtainable.  From  a  larger 
perspective,  we  have  successfully  demonstrated  the 
utility  of  a  structural  approach  to  error  detection  for 
highly  sampled  image  data.  Moreover,  we  have  built 
in  EDS  a  structure  upon  which  could  be  based  even 
more  powerful  error  detection  algorithms. 
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Appendix  A 
ALGORITHMS 


Nine  algorithms  are  provided  in  this  Appendix  (see 
Table  A-l),  in  a  block-structured  format  adopted  from 
Ref.  13. 

Individual  steps  within  each  algorithm  are  refer¬ 
enced  with  statements  of  the  form: 


Step  (X.YY.Z), 


where 


X  =  algorithm  number 
YY  =  block  number 
Z  =  step  number  within  block  YY. 

For  example,  step  2.13.1  corresponds  to  the  statement 
(from  PRMTVS): 

Set  <7,  =  ((L/LL)  +  (R/RR)  J . 

The  version  of  EDS  evaluated  statistically  in  the  text 
did  not  include  algorithms  PK2MOD  and  MERGER. 
In  addition,  step  1.7  was  inactive  during  the  numeri¬ 
cal  validations. 


Table  A-1 

Number,  name,  and  a  brief  description  of  nine 
algorithms. 


Number  Name 

1  EDS 

2  PRMTVS 

3  PK2MOD 

4  BOUNDS 

5  ISOERR 

6  NTRPL8 

7  MLTPLS 

8  SMOOTH 

9  MERGER 


Description _ 

Error  detection  and 
smoothing  (executive 
routine) 

Detects  error  peak  and  edge 
primitives 

Modifies  definition  of  PK2 
primitive 

Develops  estimates  for 
artifact  boundaries 
Detects  and  smoothes 
isolated  errors 
Single  sample  interpolation 
Detects  and  smoothes 
multiple-error  artifacts 
Linear  interpolation 
through  multi-error  artifact 
Merges  proximate  artifacts 


Algorithm  1:  EDS  Error  detection  and  smoothing,  executive  routine 

Input:  Data  points  (x(k),  k»1 ,2, . . . ,i|00) 

Output:  Smoothed  data  points  (x(k),  k= 3, A, . . . ,398} 

Note:  This  routine  uses  Algorithms  PRMTVS,  BOUNDS,  ISOERR,  and  MLTPLS. 

1.  Execute  PRMTVS. 

2.  If  all  primitives  are  null  then  terminate. 

3.  Execute  BOUNDS. 

A.  Execute  ISOERR:  Corrects  Isolated  errors. 

5.  If  any  Isolated  errors  were  found  in  step  A  then  do  block  51. 

Begin  block  5' ; 

1.  Execute  PRMTVS. 

2.  If  all  primitives  are  null  then  terminate. 

3.  Execute  BOUNDS. 

End  block  SI ; 

6.  Execute  MLTPLS:  Smoothes  multiple-error  artifarts. 

7.  If  corrected  value  differs  rrom  original  value  by  A  less, 
restore  the  pre-correction  value  1 optional). 


"T.  Pavlidis,  Structural  Pattern  Recognition,  Springer-Verlag, 
Berlin,  Heidelberg,  New  York,  1977. 
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A! go" i  trim  2  :  PRMTt'S  Detects  err o-  peak  and  edge  primitives 

Input:  Data  points  <xfk),  k»1,2 . 400} 

Output:  Primitives  { p  (k  > ,  k=J,4 . 3‘M  ;  (of.  Eq.  13) 

Note:  This  routine  uses  Algnritnn  PK2MOD. 

1.  For  k-3,A . 398  do  block  11. 

Begin  block  If  :  Identify  error  peaks; 

1.  Obtain  slopes  LL,  t.  R ,  and  RR  (cf.  Eq.  19  and  Fig.  20). 

2.  Set  p(k)-0. 

3a.  If  |L|  and  J R |  are  both  1  A,  and 
b.  if  L*R  <  0 

then  do  block  12. 

Begin  bioCK.  12; 

1.  If  LL-L  £  0,  or  if  RR-R  £  0,  or  both, 
then  p(k)  -  L/ | L j . 
else  do  block  13- 
Begin  block  13; 

1.  Set  q,  =  KL/l.L)  ♦  (R/RR )  i . 

2.  If  q,  >  Q,  -  1.85,  then  p(k)  -  L/|u|. 

3.  Execute  PK2M0D. 

End  block  13; 

End  block  12; 

End  block  1 1 ; 

2.  For  k-2,3 . 398,  do  block  21. 

Begin  block  21  :  Identify  error  edges; 

1.  Obtain  slopes  L,  C,  and  R  (cf.  Eq.  20  and  Fig.  21). 

2a.  rf  |C|  i  5,  and 
b.  If  |p(k)  j  and  | p (k* 1  )|  »  1 , 
then  do  block  22. 

Begin  block  22; 

1.  IT  L-CSO,  and  If  R/CS1 , 

then  p(k)  -  p(k*1)  »  2*C/ |C | , 
else  do  block  23. 

Begin  block  23; 

1 .  If  R-CSO  and  i f  L/CS1 , 

Then  p(k)  -  p (k* 1 }  -  2-C/|c| , 
else  do  block  24. 

Begin  block  24; 

1.  If  L,  C,  and  R  have  the  same  sign, 
then  do  block  25. 

Begin  block  25; 

1.  Set  q2  -  C/(L*R>. 

2.  If  q,  2  Q2  ■  2, 

then  p(k)  -  p(k»1)  •  2*C/|C|. 
End  block  25; 

End  block  24; 

End  block  23; 

End  block  22; 

End  block  21 ; 
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Algorithm  3:  PK2M0D  Modifies  definition  of  PK2  primitive 

Input:  a)  Data  points  !x(i),  i«1 ,2, , . . ,400i 

b)  Sample  number  K 

c)  Initial  v.alue  of  p(xl 

Output:  Revised  value  of  p<k) 

Note:  This  routine  is  used  by  Algorithm  PRMTVS. 

1.  If  | p ( k  )  |  =  1,  k  =  3.  or  k-398,  then  return. 

2.  Set  h,  *  x ( k )  -  [x(k-2)  ♦  x(k*2)]/2. 

3.  Do  block  31 . 

Begin  block  31:  assure  that  peak  is  out  of  noise; 

1.  If  h.  2  15  do  block  32; 

Begin  block  32:  extend  PK1  test  to  7-point  neighborhood; 

1.  Set  h2  -  x(k )  -  [x<k-3)  ♦  (k*3)]/2. 

2.  Set  L  =  [x(k)  -  x(k-1 ) ]. 

3.  Set  q3  =  h,/h,. 

4.  If  q,  >  Q3  -  0.8, 
then  set  p(k)  =  L/|l|. 

End  block  32; 

End  block  31 ; 
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Algorithm  A;  BOUNDS  Calculates  estimates  Tor  artifact  boundaries 
Input:  Primitives  ip(k),  k  =  3.4....,  3‘>S  } 

Output:  Estimates  for 

a)  number  of  artifacts,  -I 

b)  beginning  sample  number  of  jtn  artifact,  <b(j),  j-»1,2 . JJ 

c)  final  sample  number  of  jth  artifact,  |f!j),  j-1,2 . Ji 

1.  Set  n=0  to  Initialize  the  artifact  counter. 

2.  Set  k  =  2  to  initialize  the  p(k)  pointer. 

3.  While  k<39S  do  block  11. 


Begin  block  1 1 ; 

1.  Set  k  -  k*1  to  move  the  p(k)  pointer  along  the  scan  line. 

2.  If  p(k )  «  0  do  block  12. 


Begin  block  12:  an  artifact  has  been  found; 

1.  Set  n  -  n»1  to  increment  the  artifact  counter. 

2.  Set  b(n)  =  k  and  FLAG  -  1. 

3.  Do  block  13  for  i  -  k,  395. 

Begin  block  13:  search  for  end  of  artifact; 

1 .  Set  f (n)  -  1. 

2.  If  p[f (n)*1 ]  =  p[f (n)+2]  =  p[f(n)+31  -  0  go  to  12.7 
End  block  13; 

A.  Set  f (n)  -  398. 

5.  If  p(398)  >  0  set  f(n)  -  397. 

6.  If  p(397)  -  p(398)  -  0  set  f(n)  -  396. 

7.  If  b(n )  S  A  or  f(n)  i  397  set  FLAG  -  0. 

8.  ir  p[b(n) ]  •  p[r(n)]  •  FLAG  -  -A  do  block  1A. 

Begin  block  1A:  shrink  plateau  artifact  boundaries; 

1 .  Set  b(n)  »  b(n)  ♦  1 . 

2.  Set  f (n)  -  f (n)  -  1 . 

End  block  1A; 

9.  Set  k  =  f(n)+3;  search  for  next  artifact  will  start  at 
f(n)»A. 

End  block  12; 

End  block  1 1 ; 
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Algorithm  5:  ISOF.RR  Detect:;  jml  smoothes  isolated  errors 

Input:  a)  Data  points  (x(k),  <  =  1,2 . *100) 

b)  Artifact  boundaries  lb(jl,  f(j',  j*1 .2 . J ) 

Output:  Data  points  {x'k',  k«l,2 . 400i  ,  with  isolated  errors 

corrected  over  the  range  3  S  k  <  398. 

Note:  This  routine  uses  Algorithms  NTRPL8  and  PRMTVS. 

1 .  Do  block  1 1  for  j-1 ,2. 

Begin  block  11; 

1.  Set  UJ )  »  [f(j)  -  b(j)  »  1  ]  to  estimate  artifact  length. 

2.  If  L(j)  S  4  do  block  12. 

Begin  block  12  :  Attempt  single  sample  correction  of  artifact  J 

1.  Set  i=0,  B(j)  =  b( j ) ,  and  F(j)  -  f 'j ). 

2.  If  L(J)-1  do  block  125. 

Begin  block  125  :  Expand  search  region; 

1.  If  B(j)  >  3  set  B(  j  1  =  [  B(  j  )  -  1]. 

2.  If  F  ( j  )  <  398  set  FCJ )  -  [FCj)  ♦ 

End  block  125; 

3.  Do  block  13  for  n  =  B( j  > ,  Fiji. 

Begin  block  13  :  Attempt  single  sample  correction  at  n; 

1.  Set  2  =*  x(n)  to  save  original  value. 

2.  Apply  NTRPL8  to  (x(k)J  to  get  y(n). 

3.  Set  x(n)  -  y(n)  to  perform  replacement. 

4.  Apply  PRMTVS  to  (x(k)}  to  get  ip’(k)}. 

5.  If  p'(k).0,  BSkSF,  do  block  135. 

Begin  block  135  :  Single  correction  at  n  is 
viable ; 

1 .  Set  i-i  +  1  . 

2.  Set  Nil)  -  n. 

3.  Calculate  a(i)  from  Sq.  >22). 

End  block  135; 

6.  Set  x(n)  »  z  to  restore  original  value. 

End  block  13; 

4.  If  1.0  do  block  14. 

Begin  block  14  :  Apply  best  of  the  i  viable  corrections. 

1.  Find  m  such  that  a(m)  -  min{a(1),  a(2) . a(i)). 

2.  Set  x[N(m)])  *  y[N(m)],  to  perform  correction. 

End  block  14; 

End  block  12; 

End  block  1 1  : 
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De f  i -i i  tlO'„s  of  Vj-  (abl""  in  A : f.if  1  tnn  T.'-O-'RR 


j  indexes  tit'  "current  i'*tifart",  i.e.,  the  u-'tifact  cur'-enl'y 

being  considered  for  single  sample  correct : on. 

J  maxi-nun  v  t  i  a-  of  j;  octal  number  of  jrti  facts  in  toe  waveform 

(input  to  Algo'*;  t.nm  ISOSRR .  provide  !  by  Algoritnm  sOUSDS) . 

9(j),F(j)  bound  tie  argument  range  of  <!x)  over  which  the  search  will  oe 
conducted  for  a  single  sample  correction  to  artifact  j. 

n  running  index  of  position  within  cu*’rent  artifact. 

z  dummy  variable  for  temporarily  storing  a  pre-correction  value 

of  x ( x  ) . 

y(n)  a  -eplacement  value  for  x(n),  generated  by  interpolation  from 

adjacent  samples  (output  of  Algorithm  NTRPL8). 

p • ( k  )  sequence  of  error  primitives  for  waveform  In  which  sample  x(n!  has 

been  replaced  by  interpolated  value  (output  of  Algorithm  PRMTVS). 

i  counts  the  number  of  viable  single  sample  corrections  within  the 

current  artifact. 

N(i)  location  of  the  ith  viable  single  sample  correction  within  the 

current  artifact. 

a(i)  measure  of  smoothness  of  sequence  !x(k)i  after  ith  viable 

Isolated  correction  to  the  current  artifact  (of.  Eq.  22). 

m  value  of  i  tnat  indexes  the  "best"  of  the  viable  isolated 

corrections,  in  that  the  resulting  corrected  sequence  { x ( < ) } 
is  "smoothest". 


Algorithm  6;  NTRPL.8  Single  sample  Interpolation 

Input:  a)  Data  points  (x(k),  k=1 ,2, . . . ,400} 

b)  Location,  n,  of  sample  for  which  replacement  is  desired. 

Output:  Interpolated  value,  y(n>,  to  replace  the  presumably  defective 

value,  x(n). 

Note:  This  routine  is  used  by  Algorithm  ISOERR. 

Integer  n  is  always  within  the  range  3SnS398. 


1.  Set  y(n)  -  [x(n-1)  ♦  x(n»l)]/2,  to  obtain  linear  interpolation. 

2.  Set  C  -  [x(n-2)  *  x(n*2)]/2,  to  obtain  correction  term. 

3.  Set  y(n)  »  y(n)  ♦  [y(n)-C]/3.  to  obtain  cubic  interpolation. 
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Algorithm  7:  MLTPLS  Detects  and  smoothes  multiple-error  artifacts 

Input:  a)  Data  points  !x(k),  k-1 ADO!  with  isolated  errors  already 

corrected. 

b)  Primitives  { p  ( k ) ,  k=3. 4 . 398} 

c)  Estimated  number  of  multiple-error  artifacts  J 

d)  Estimated  artifact  boundaries  {b(j>,  f(j),  j«1,?,...,J) 

Output:  Smoothed  data  points  !x(kl,  k»3«4, . . . , 398) 

Mote:  This  algorithm  uses  Algorithms  PRMTi/S .  SMOOTH,  and  MERGER. 


J.  Set  j-0  and  n=  0. 

2.  While  j  <  J  do  block  1 1 . 

Begin  block  11; 

1.  Set  j=J  +  1  and  n»n*-1. 

2.  Set  B<n)  -  S(j )  and  F(n)  -  f(j). 

3.  Set  m«1 . 

«.  Set  p'(2)  -  p'(399)-0. 

5.  If  m*0  do  block  12. 

Begin  block  12; 

1 .  Set  m»0. 

2.  Execute  SMOOTH  to  interpolate  linearly  through  multi-error 
artifact. 

3-  Apply  PRMTVS  to  (x(k)l  to  obtain  (p'(k)h 

4.  Do  block  21  for  k  «  [B(n)-1],  [F(n)+1]. 

Begin  block  21:  Count  the  non-null  entries  in  p'(k); 

1.  If  p'M*0  do  block  22. 

Begin  block  22:  Error  artifact  still  exists; 

1 .  Set  m*m*1 . 

2.  Set  km-k. 

End  block  22; 

End  block  21 ; 

5.  If  m*0  do  block  31. 

Begin  block  31:  Artifact  still  exl3t3  -  grow  the 
boundaries: 

1.  Set  B(n)  -  minfB(n) ,k, ) . 

2.  Set  F(n)  -  max  (F(n),km). 

3.  Execute  MERGER  to  merge  proximate  artifacts. 

End  block  31 ; 

End  block  12; 

End  block  1 1 ; 


THE  JOHNS  HOPKINS  UNIVERSITY 

APPLIED  PHYSICS  LABORATORY 

LAUREL.  MARYLAND 


-TV 


Algorithm  3:  SMOOTH  Linear  interpolation  thr.ig'i  m  il  t; -error  d-tifaet 

Input:  a)  Index  n 

b)  Data  points  [*!'<),  <» I ,?,... ,k00) 

c)  Boundaries  (3(i),  F(i',  1*1,2 . n) 


Output : 


Smoothed  data  { x ( !<  > ,  x-1,2. 


This  routine  is  used  by  Algo-ir.hm  MLIPLS. 


1.  If  Bin  )S5  do  block  11  for  k  *1 , . . .  ,F(n> . 

Begin  block  11:  blank  beginning  of  line; 

1 .  Set  x(k)  -  x[F(n)*1 ]. 

End  block  1 1 ; 

2.  If  F(n)i396  do  block  21  for  k*B(n1 . U00. 

Begin  block  21:  blank  end  of  line; 

1.  Set  x(k)  *  x[B(n)-1]. 

End  block  21 ; 

3.  If  B(n)s6  and  F(n)S395  do  block  31. 

Begin  block  31:  linear  interpolation  from  (B-1)  to  (F*1); 

1.  Set  L(n)  -  F(n )  -  B(n)  ♦  1. 

2.  Set  d  -  (x[F(.n)*1]  -  x[B(n)-1  ]}/'L(n)*1  ]. 

3.  Do  block  32  for  i-1,  L(n). 

Begin  block  32 j 

1.  Set  k  -  B(n)*i-1. 

2.  Set  x(k)  -  x[B(n)-1 ]  ♦  i*S. 

End  block  32; 

End  block  31 ; 


Algorithm  9:  MERGER  Merges  proximate  artifacts 

Input:  a)  Indices  J  and  n 

b)  Boundary  estimates  lb(i),  f(i),  i*1,2,...J) 

c)  Improved  boundaries  (B(t),  F(i),  1-1,2, ...n! 

Output:  a)  Revised  values  of  j  and  n 

b)  Revised  values  of  Bin)  and  F(n) 

Note:  This  routine  is  used  by  Algorithm  MLTPLS. 

1 .  Do  block  1 1 . 

Begin  block  1 1 ; 

1 .  Set  C.-2  and  Cg*5. 

2.  Set  l(n)  -  F(n)-B(n)*l. 

3.  If  L(n)iC2  set  Cj-CjM. 

End  block  11; 

2.  If  JS(J-I)  and  F(n)i[b(J*1 )  -  C,]  then  do  block  21. 

Begin  block  21:  Region  merging,  forward; 

1.  Set  J-J*1. 

2.  Set  F(n)  -  f(J). 

End  block  21 ; 

3.  If  ni2  and  B(n)S[F(n)  ♦  C,]  then  do  block  31. 

Begin  block  31:  Region  merging,  backward; 

1 .  Set  n-n-1 . 

2.  Set  F(n )  -  F(n*1). 

End  block  31 ; 


vm 
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DERIVATION  OF  CUBIC  INTERPOLATION  FORMULA 


Algorithm  NTRPL8  generates  a  replacement  value 
for  x(n),  called y(n),  by  interpolating  from  the  four 
neighboring  values:  jr(/>  -  2),x(n  -  \),x (/»  +  1), 
and  x(n  +  2).  In  this  appendix  we  derive  the  interpo¬ 
lation  formula  used  by  NTRPL8. 

Since  it  is  our  intent  to  interpolate  between  four 
“known”  data  points,  we  select  as  our  interpolation 
function  a  cubic  polynomial, 

y(k)  =  a  +  b(k  ~  n)  +  c(k  -  n)2 

+  d(k  -  n)\  Ik  -  n\  £  2,  (B-l) 

where  the  coefficients  (a,  b,  c,  and  d)  are  at  present 
treated  as  “unknowns”  in  need  of  determination. 

The  four  data  values  neighboring  sample  n  are  sub¬ 
stituted  into  Eq.  B-l  to  obtain  four  equations: 


x(n  -  2)  =  y(n  -  2) 

=  o-  2b  +  4c-8d  (B-2) 

x(.n  -  I)  =  y(n  -  1) 

=  a  —  b  +  c  —  d  (B-3) 

x(n  +  1)  =  y(n  +  1) 

=  a  +  b  +  c  +  d  (B-4) 

x(n  +  2)  =  y(n  +  2) 

=  a  +  2b +  4c  +  Sd.  (B-5) 


We  next  define  quantities  St  and  C  as  follows: 
k  -  lx(n  -  1)  +  x(n  +  1)1/2,  (B-6) 


From  Eqs.  B-3,  B-4,  and  B-6, 

St  -  a  +  c  .  (B-8) 

From  Eqs.  B-2,  B-S,  and  B-7, 

C  =  a  +  4c  .  (B-9) 

From  Eqs.  B-8  and  B-9, 

a  =  (4£  -  C)/3  =  St  +  (*  -  0/3  .  (B-10) 

However,  from  Eq.  B-l, 

y(n)  -  a  .  (B-l  1) 

Thus,  from  Eqs.  B-10  and  B-ll, 

y(n)  *  St  +  (St  -  0/3  .  (B-I2) 

where 

St  m  (jr(n  -  1)  +  x(n  +  1)1/2  (B-13) 

C  =  [jr(n  -  2)  +  x(n  +  2)1/2  .  (B-14) 

The  linear  interpolation  estimate,  Eq.  B-13,  is 
formed  in  step  6.1;  the  correction  term  C,  Eq.  B-14, 
is  formed  in  step  6.2;  the  cubic  interpolation,  Eq.  B-12, 
is  formed  in  step  6.3. 

It  also  follows  from  Eqs.  B-12  through  B-14  that 


& 


and 


C  =  lx  (ft  -  2)  +  x(n  +  2)1/2  .  (B-7) 


y(n)  -  (-•*■(«  -  2)  +  4x(n  -  1)  + 

4x(n  +  1  )-*(/»  +  2)J/6  .  (B-15) 
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